1 | MODULE domzgr |
---|
2 | !!============================================================================== |
---|
3 | !! *** MODULE domzgr *** |
---|
4 | !! Ocean initialization : domain initialization |
---|
5 | !!============================================================================== |
---|
6 | |
---|
7 | !!---------------------------------------------------------------------- |
---|
8 | !! dom_zgr : defined the ocean vertical coordinate system |
---|
9 | !! zgr_bat : bathymetry fields (levels and meters) |
---|
10 | !! zgr_bat_zoom : modify the bathymetry field if zoom domain |
---|
11 | !! zgr_bat_ctl : check the bathymetry files |
---|
12 | !! zgr_z : reference z-coordinate |
---|
13 | !! zgr_zps : z-coordinate with partial steps |
---|
14 | !! zgr_s : s-coordinate |
---|
15 | !!--------------------------------------------------------------------- |
---|
16 | !! * Modules used |
---|
17 | USE oce ! ocean dynamics and tracers |
---|
18 | USE dom_oce ! ocean space and time domain |
---|
19 | USE in_out_manager ! I/O manager |
---|
20 | USE lib_mpp ! distributed memory computing library |
---|
21 | USE lbclnk ! ocean lateral boundary conditions (or mpp link) |
---|
22 | USE closea |
---|
23 | USE solisl |
---|
24 | |
---|
25 | IMPLICIT NONE |
---|
26 | PRIVATE |
---|
27 | |
---|
28 | !! * Routine accessibility |
---|
29 | PUBLIC dom_zgr ! called by dom_init.F90 |
---|
30 | |
---|
31 | !! * Substitutions |
---|
32 | # include "domzgr_substitute.h90" |
---|
33 | # include "vectopt_loop_substitute.h90" |
---|
34 | !!---------------------------------------------------------------------- |
---|
35 | !! OPA 9.0 , LODYC-IPSL (2003) |
---|
36 | !!---------------------------------------------------------------------- |
---|
37 | |
---|
38 | CONTAINS |
---|
39 | |
---|
40 | SUBROUTINE dom_zgr |
---|
41 | !!---------------------------------------------------------------------- |
---|
42 | !! *** ROUTINE dom_zgr *** |
---|
43 | !! |
---|
44 | !! ** Purpose : set the depth of model levels and the resulting |
---|
45 | !! vertical scale factors. |
---|
46 | !! |
---|
47 | !! ** Method reference vertical coordinate |
---|
48 | !! Z-coordinates : The depth of model levels is defined |
---|
49 | !! from an analytical function the derivative of which gives |
---|
50 | !! the vertical scale factors. |
---|
51 | !! both depth and scale factors only depend on k (1d arrays). |
---|
52 | !! w-level: gdepw = fsdep(k) |
---|
53 | !! e3w(k) = dk(fsdep)(k) = fse3(k) |
---|
54 | !! t-level: gdept = fsdep(k+0.5) |
---|
55 | !! e3t(k) = dk(fsdep)(k+0.5) = fse3(k+0.5) |
---|
56 | !! |
---|
57 | !! ** Action : - gdept, gdepw : depth of T- and W-point (m) |
---|
58 | !! - e3t, e3w : scale factors at T- and W-levels (m) |
---|
59 | !! |
---|
60 | !! Reference : |
---|
61 | !! Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766. |
---|
62 | !! |
---|
63 | !! History : |
---|
64 | !! 9.0 ! 03-08 (G. Madec) original code |
---|
65 | !!---------------------------------------------------------------------- |
---|
66 | INTEGER :: ioptio = 0 ! temporary integer |
---|
67 | !!---------------------------------------------------------------------- |
---|
68 | |
---|
69 | ! Check Vertical coordinate options |
---|
70 | ! --------------------------------- |
---|
71 | ioptio = 0 |
---|
72 | IF( lk_sco ) THEN |
---|
73 | IF(lwp) WRITE(numout,*) |
---|
74 | IF(lwp) WRITE(numout,*) 'dom_zgr : s-coordinate' |
---|
75 | IF(lwp) WRITE(numout,*) '~~~~~~~' |
---|
76 | ioptio = ioptio + 1 |
---|
77 | ENDIF |
---|
78 | IF( lk_zps ) THEN |
---|
79 | IF(lwp) WRITE(numout,*) |
---|
80 | IF(lwp) WRITE(numout,*) 'dom_zgr : z-coordinate with partial steps' |
---|
81 | IF(lwp) WRITE(numout,*) '~~~~~~~' |
---|
82 | ioptio = ioptio + 1 |
---|
83 | ENDIF |
---|
84 | IF( ioptio == 0 ) THEN |
---|
85 | IF(lwp) WRITE(numout,*) |
---|
86 | IF(lwp) WRITE(numout,*) 'dom_zgr : z-coordinate' |
---|
87 | IF(lwp) WRITE(numout,*) '~~~~~~~' |
---|
88 | ENDIF |
---|
89 | |
---|
90 | IF ( ioptio > 1 ) THEN |
---|
91 | IF(lwp) WRITE(numout,cform_err) |
---|
92 | IF(lwp) WRITE(numout,*) ' several vertical coordinate options used' |
---|
93 | nstop = nstop + 1 |
---|
94 | ENDIF |
---|
95 | |
---|
96 | ! Build the vertical coordinate system |
---|
97 | ! ------------------------------------ |
---|
98 | |
---|
99 | CALL zgr_z ! Reference z-coordinate system |
---|
100 | |
---|
101 | CALL zgr_bat ! Bathymetry fields (levels and meters) |
---|
102 | |
---|
103 | CALL zgr_zps ! Partial step z-coordinate |
---|
104 | |
---|
105 | CALL zgr_s ! s-coordinate |
---|
106 | |
---|
107 | END SUBROUTINE dom_zgr |
---|
108 | |
---|
109 | |
---|
110 | SUBROUTINE zgr_z |
---|
111 | !!---------------------------------------------------------------------- |
---|
112 | !! *** ROUTINE zgr_z *** |
---|
113 | !! |
---|
114 | !! ** Purpose : set the depth of model levels and the resulting |
---|
115 | !! vertical scale factors. |
---|
116 | !! |
---|
117 | !! ** Method : z-coordinate system (use in all type of coordinate) |
---|
118 | !! The depth of model levels is defined from an analytical |
---|
119 | !! function the derivative of which gives the scale factors. |
---|
120 | !! both depth and scale factors only depend on k (1d arrays). |
---|
121 | !! w-level: gdepw = fsdep(k) |
---|
122 | !! e3w(k) = dk(fsdep)(k) = fse3(k) |
---|
123 | !! t-level: gdept = fsdep(k+0.5) |
---|
124 | !! e3t(k) = dk(fsdep)(k+0.5) = fse3(k+0.5) |
---|
125 | !! |
---|
126 | !! ** Action : - gdept, gdepw : depth of T- and W-point (m) |
---|
127 | !! - e3t, e3w : scale factors at T- and W-levels (m) |
---|
128 | !! |
---|
129 | !! Reference : |
---|
130 | !! Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766. |
---|
131 | !! |
---|
132 | !! History : |
---|
133 | !! 9.0 ! 03-08 (G. Madec) F90: Free form and module |
---|
134 | !!---------------------------------------------------------------------- |
---|
135 | !! * Local declarations |
---|
136 | INTEGER :: jk ! dummy loop indices |
---|
137 | REAL(wp) :: zt, zw ! temporary scalars |
---|
138 | REAL(wp) :: & |
---|
139 | zsur , za0, za1, zkth, zacr, & ! Values set from parameters in |
---|
140 | zdzmin, zhmax ! par_CONFIG_Rxx.h90 |
---|
141 | !!---------------------------------------------------------------------- |
---|
142 | |
---|
143 | ! Set variables from parameters |
---|
144 | ! ------------------------------ |
---|
145 | zkth = ppkth ; zacr = ppacr |
---|
146 | zdzmin = ppdzmin ; zhmax = pphmax |
---|
147 | |
---|
148 | ! If ppa1 and ppa0 and ppsur are et to pp_to_be_computed |
---|
149 | ! za0, za1, zsur are computed from ppdzmin , pphmax, ppkth, ppacr |
---|
150 | ! |
---|
151 | IF( ppa1 == pp_to_be_computed .AND. & |
---|
152 | & ppa0 == pp_to_be_computed .AND. & |
---|
153 | & ppsur == pp_to_be_computed ) THEN |
---|
154 | za1 = ( ppdzmin - pphmax / FLOAT(jpk-1) ) & |
---|
155 | / ( TANH((1-ppkth)/ppacr) - ppacr/FLOAT(jpk-1) & |
---|
156 | & * ( LOG( COSH( (jpk - ppkth) / ppacr) ) & |
---|
157 | & - LOG( COSH( ( 1 - ppkth) / ppacr) ) ) ) |
---|
158 | |
---|
159 | za0 = ppdzmin - za1 * TANH( (1-ppkth) / ppacr ) |
---|
160 | zsur = - za0 - za1 * ppacr * LOG( COSH( (1-ppkth) / ppacr ) ) |
---|
161 | |
---|
162 | ELSE |
---|
163 | za1 = ppa1 ; za0 = ppa0 ; zsur = ppsur |
---|
164 | ENDIF |
---|
165 | |
---|
166 | |
---|
167 | ! Parameter print |
---|
168 | ! --------------- |
---|
169 | IF(lwp) THEN |
---|
170 | WRITE(numout,*) |
---|
171 | WRITE(numout,*) ' zgr_z : Reference vertical z-coordinates' |
---|
172 | WRITE(numout,*) ' ~~~~~~~' |
---|
173 | IF ( ppa1 == 0. .AND. ppa0 == 0. .AND. ppsur == 0. ) THEN |
---|
174 | WRITE(numout,*) ' zsur, za0, za1 computed from ' |
---|
175 | WRITE(numout,*) ' zdzmin = ', zdzmin |
---|
176 | WRITE(numout,*) ' zhmax = ', zhmax |
---|
177 | ENDIF |
---|
178 | WRITE(numout,*) ' Namelist namzgr : value of coefficients for vertical mesh:' |
---|
179 | WRITE(numout,*) ' zsur = ', zsur |
---|
180 | WRITE(numout,*) ' za0 = ', za0 |
---|
181 | WRITE(numout,*) ' za1 = ', za1 |
---|
182 | WRITE(numout,*) ' zkth = ', zkth |
---|
183 | WRITE(numout,*) ' zacr = ', zacr |
---|
184 | ENDIF |
---|
185 | |
---|
186 | |
---|
187 | ! Reference z-coordinate (depth - scale factor at T- and W-points) |
---|
188 | ! ====================== |
---|
189 | |
---|
190 | DO jk = 1, jpk |
---|
191 | zw = FLOAT( jk ) |
---|
192 | zt = FLOAT( jk ) + 0.5 |
---|
193 | gdepw(jk) = ( zsur + za0 * zw + za1 * zacr * LOG( COSH( (zw-zkth)/zacr ) ) ) |
---|
194 | gdept(jk) = ( zsur + za0 * zt + za1 * zacr * LOG( COSH( (zt-zkth)/zacr ) ) ) |
---|
195 | e3w (jk) = za0 + za1 * TANH( (zw-zkth)/zacr ) |
---|
196 | e3t (jk) = za0 + za1 * TANH( (zt-zkth)/zacr ) |
---|
197 | END DO |
---|
198 | gdepw(1) = 0.e0 ! force first w-level to be exactly at zero |
---|
199 | |
---|
200 | |
---|
201 | ! Control and print |
---|
202 | ! ================== |
---|
203 | |
---|
204 | IF(lwp) THEN |
---|
205 | WRITE(numout,*) |
---|
206 | WRITE(numout,*) ' Reference z-coordinate depth and scale factors:' |
---|
207 | WRITE(numout, "(9x,' level gdept gdepw e3t e3w ')" ) |
---|
208 | WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, gdept(jk), gdepw(jk), e3t(jk), e3w(jk), jk = 1, jpk ) |
---|
209 | ENDIF |
---|
210 | |
---|
211 | DO jk = 1, jpk |
---|
212 | IF( e3w(jk) <= 0. .OR. e3t(jk) <= 0. ) THEN |
---|
213 | IF(lwp) WRITE(numout,cform_err) |
---|
214 | IF(lwp) WRITE(numout,*) ' e3w or e3t =< 0 ' |
---|
215 | nstop = nstop + 1 |
---|
216 | ENDIF |
---|
217 | IF( gdepw(jk) < 0. .OR. gdept(jk) < 0.) THEN |
---|
218 | IF(lwp) WRITE(numout,cform_err) |
---|
219 | IF(lwp) WRITE(numout,*) ' gdepw or gdept < 0 ' |
---|
220 | nstop = nstop + 1 |
---|
221 | ENDIF |
---|
222 | END DO |
---|
223 | |
---|
224 | END SUBROUTINE zgr_z |
---|
225 | |
---|
226 | |
---|
227 | SUBROUTINE zgr_bat |
---|
228 | !!---------------------------------------------------------------------- |
---|
229 | !! *** ROUTINE zgr_bat *** |
---|
230 | !! |
---|
231 | !! ** Purpose : set bathymetry both in levels and meters |
---|
232 | !! |
---|
233 | !! ** Method : read or define mbathy and bathy arrays |
---|
234 | !! * level bathymetry: |
---|
235 | !! The ocean basin geometry is given by a two-dimensional array, |
---|
236 | !! mbathy, which is defined as follow : |
---|
237 | !! mbathy(ji,jj) = 1, ..., jpk-1, the number of ocean level |
---|
238 | !! at t-point (ji,jj). |
---|
239 | !! = 0 over the continental t-point. |
---|
240 | !! = -n over the nth island t-point. |
---|
241 | !! The array mbathy is checked to verified its consistency with |
---|
242 | !! model option. in particular: |
---|
243 | !! mbathy must have at least 1 land grid-points (mbathy<=0) |
---|
244 | !! along closed boundary. |
---|
245 | !! mbathy must be cyclic IF jperio=1. |
---|
246 | !! mbathy must be lower or equal to jpk-1. |
---|
247 | !! isolated ocean grid points are suppressed from mbathy |
---|
248 | !! since they are only connected to remaining |
---|
249 | !! ocean through vertical diffusion. |
---|
250 | !! ntopo=-1 : rectangular channel or bassin with a bump |
---|
251 | !! ntopo= 0 : flat rectangular channel or basin |
---|
252 | !! ntopo= 1 : mbathy is read in 'bathy_level.nc' NetCDF file |
---|
253 | !! bathy is read in 'bathy_meter.nc' NetCDF file |
---|
254 | !! C A U T I O N : mbathy will be modified during the initializa- |
---|
255 | !! tion phase to become the number of non-zero w-levels of a water |
---|
256 | !! column, with a minimum value of 1. |
---|
257 | !! |
---|
258 | !! ** Action : - mbathy: level bathymetry (in level index) |
---|
259 | !! - bathy : meter bathymetry (in meters) |
---|
260 | !! |
---|
261 | !! History : |
---|
262 | !! 9.0 ! 03-08 (G. Madec) Original code |
---|
263 | !!---------------------------------------------------------------------- |
---|
264 | !! * Modules used |
---|
265 | USE ioipsl |
---|
266 | |
---|
267 | !! * Local declarations |
---|
268 | CHARACTER (len=15) :: clname ! temporary characters |
---|
269 | LOGICAL :: llbon ! check the existence of bathy files |
---|
270 | INTEGER :: ji, jj, jl, jk ! dummy loop indices |
---|
271 | INTEGER :: inum = 11 ! temporary logical unit |
---|
272 | INTEGER :: & |
---|
273 | ipi, ipj, ipk, & ! temporary integers |
---|
274 | itime, & ! " " |
---|
275 | ii_bump, ij_bump ! bump center position |
---|
276 | INTEGER, DIMENSION (1) :: istep |
---|
277 | INTEGER , DIMENSION(jpidta,jpjdta) :: & |
---|
278 | idta ! global domain integer data |
---|
279 | REAL(wp) :: & |
---|
280 | r_bump, h_bump, h_oce, & ! bump characteristics |
---|
281 | zi, zj, zdate0, zdt ! temporary scalars |
---|
282 | REAL(wp), DIMENSION(jpidta,jpjdta) :: & |
---|
283 | zlamt, zphit, & ! temporary workspace (NetCDF read) |
---|
284 | zdta ! global domain scalar data |
---|
285 | REAL(wp), DIMENSION(jpk) :: & |
---|
286 | zdept ! temporary workspace (NetCDF read) |
---|
287 | !!---------------------------------------------------------------------- |
---|
288 | |
---|
289 | IF(lwp) WRITE(numout,*) |
---|
290 | IF(lwp) WRITE(numout,*) ' zgr_bat : defines level and meter bathymetry' |
---|
291 | IF(lwp) WRITE(numout,*) ' ~~~~~~~' |
---|
292 | |
---|
293 | ! ======================================== |
---|
294 | ! global domain level and meter bathymetry (idta,zdta) |
---|
295 | ! ======================================== |
---|
296 | ! ! =============== ! |
---|
297 | IF( ntopo == 0 .OR. ntopo == -1 ) THEN ! defined by hand ! |
---|
298 | ! ! =============== ! |
---|
299 | |
---|
300 | IF( ntopo == 0 ) THEN ! flat basin |
---|
301 | |
---|
302 | IF(lwp) WRITE(numout,*) |
---|
303 | IF(lwp) WRITE(numout,*) ' bathymetry field: flat basin' |
---|
304 | |
---|
305 | idta(:,:) = jpkm1 ! flat basin |
---|
306 | zdta(:,:) = gdepw(jpk) |
---|
307 | |
---|
308 | ELSE ! bump |
---|
309 | IF(lwp) WRITE(numout,*) |
---|
310 | IF(lwp) WRITE(numout,*) ' bathymetry field: flat basin with a bump' |
---|
311 | |
---|
312 | ii_bump = jpidta / 3 + 3 ! i-index of the bump center |
---|
313 | ij_bump = jpjdta / 2 ! j-index of the bump center |
---|
314 | r_bump = 6 ! bump radius (index) |
---|
315 | h_bump = 240.e0 ! bump height (meters) |
---|
316 | h_oce = gdepw(jpk) ! background ocean depth (meters) |
---|
317 | IF(lwp) WRITE(numout,*) ' bump characteristics: ' |
---|
318 | IF(lwp) WRITE(numout,*) ' bump center (i,j) = ', ii_bump, ii_bump |
---|
319 | IF(lwp) WRITE(numout,*) ' bump height = ', h_bump , ' meters' |
---|
320 | IF(lwp) WRITE(numout,*) ' bump radius = ', r_bump , ' index' |
---|
321 | IF(lwp) WRITE(numout,*) ' background ocean depth = ', h_oce , ' meters' |
---|
322 | ! zdta : |
---|
323 | DO jj = 1, jpjdta |
---|
324 | DO ji = 1, jpidta |
---|
325 | zi = FLOAT( ji - ii_bump ) / r_bump |
---|
326 | zj = FLOAT( jj - ij_bump ) / r_bump |
---|
327 | zdta(ji,jj) = h_oce - h_bump * EXP( -( zi*zi + zj*zj ) ) |
---|
328 | END DO |
---|
329 | END DO |
---|
330 | ! idta : |
---|
331 | idta(:,:) = jpkm1 |
---|
332 | DO jk = 1, jpkm1 |
---|
333 | DO jj = 1, jpjdta |
---|
334 | DO ji = 1, jpidta |
---|
335 | IF( gdept(jk) < zdta(ji,jj) .AND. zdta(ji,jj) <= gdept(jk+1) ) idta(ji,jj) = jk |
---|
336 | END DO |
---|
337 | END DO |
---|
338 | END DO |
---|
339 | ENDIF |
---|
340 | |
---|
341 | ! set boundary conditions (caution, idta on the global domain: use of jperio, not nperio) |
---|
342 | IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) THEN |
---|
343 | idta( : , 1 ) = -1 ; zdta( : , 1 ) = -1.e0 |
---|
344 | idta( : ,jpjdta) = 0 ; zdta( : ,jpjdta) = 0.e0 |
---|
345 | ELSEIF( jperio == 2 ) THEN |
---|
346 | idta( : , 1 ) = idta( : , 3 ) ; zdta( : , 1 ) = zdta( : , 3 ) |
---|
347 | idta( : ,jpjdta) = 0 ; zdta( : ,jpjdta) = 0.e0 |
---|
348 | idta( 1 , : ) = 0 ; zdta( 1 , : ) = 0.e0 |
---|
349 | idta(jpidta, : ) = 0 ; zdta(jpidta, : ) = 0.e0 |
---|
350 | ELSE |
---|
351 | idta( : , 1 ) = 0 ; zdta( : , 1 ) = 0.e0 |
---|
352 | idta( : ,jpjdta) = 0 ; zdta( : ,jpjdta) = 0.e0 |
---|
353 | idta( 1 , : ) = 0 ; zdta( 1 , : ) = 0.e0 |
---|
354 | idta(jpidta, : ) = 0 ; zdta(jpidta, : ) = 0.e0 |
---|
355 | ENDIF |
---|
356 | |
---|
357 | ! EEL R5 configuration with east and west open boundaries. |
---|
358 | ! Two rows of zeroes are needed at the south and north for OBCs |
---|
359 | |
---|
360 | IF( cp_cfg == "eel" .AND. jp_cfg == 5 ) THEN |
---|
361 | idta( : , 2 ) = 0 ; zdta( : , 2 ) = 0.e0 |
---|
362 | idta( : ,jpjdta-1) = 0 ; zdta( : ,jpjdta-1) = 0.e0 |
---|
363 | ENDIF |
---|
364 | |
---|
365 | ! ! =============== ! |
---|
366 | ELSEIF( ntopo == 1 ) THEN ! read in file ! |
---|
367 | ! ! =============== ! |
---|
368 | |
---|
369 | clname = 'bathy_level.nc' ! Level bathymetry |
---|
370 | INQUIRE( FILE=clname, EXIST=llbon ) |
---|
371 | IF( llbon ) THEN |
---|
372 | IF(lwp) WRITE(numout,*) |
---|
373 | IF(lwp) WRITE(numout,*) ' read level bathymetry in ', clname |
---|
374 | IF(lwp) WRITE(numout,*) |
---|
375 | itime = 1 |
---|
376 | ipi = jpidta |
---|
377 | ipj = jpjdta |
---|
378 | ipk = 1 |
---|
379 | zdt = rdt |
---|
380 | CALL flinopen( clname, 1, jpidta, 1, jpjdta, .FALSE., & |
---|
381 | ipi, ipj, ipk, zlamt, zphit, zdept, itime, istep, zdate0, zdt, inum ) |
---|
382 | CALL flinget( inum, 'Bathy_level', jpidta, jpjdta, 1, & |
---|
383 | itime, 1, 1, 1, jpidta, 1, jpjdta, zdta(:,:) ) |
---|
384 | idta(:,:) = zdta(:,:) |
---|
385 | CALL flinclo( inum ) |
---|
386 | |
---|
387 | ELSE |
---|
388 | IF(lwp) WRITE(numout,cform_err) |
---|
389 | IF(lwp) WRITE(numout,*)' zgr_bat : unable to read the file', clname |
---|
390 | nstop = nstop + 1 |
---|
391 | ENDIF |
---|
392 | |
---|
393 | clname = 'bathy_meter.nc' ! meter bathymetry |
---|
394 | INQUIRE( FILE=clname, EXIST=llbon ) |
---|
395 | IF( llbon ) THEN |
---|
396 | IF(lwp) WRITE(numout,*) |
---|
397 | IF(lwp) WRITE(numout,*) ' read meter bathymetry in ', clname |
---|
398 | IF(lwp) WRITE(numout,*) |
---|
399 | itime = 1 |
---|
400 | ipi = jpidta |
---|
401 | ipj = jpjdta |
---|
402 | ipk = 1 |
---|
403 | zdt = rdt |
---|
404 | CALL flinopen( clname, 1, jpidta, 1, jpjdta, .FALSE., & |
---|
405 | ipi, ipj, ipk, zlamt, zphit, zdept, itime, istep, zdate0, zdt, inum ) |
---|
406 | CALL flinget( inum, 'Bathymetry', jpidta, jpjdta, 1, & |
---|
407 | itime, 1, 1, 1, jpidta, 1, jpjdta, zdta(:,:) ) |
---|
408 | CALL flinclo( inum ) |
---|
409 | ELSE |
---|
410 | IF( lk_zps .OR. lk_sco ) THEN |
---|
411 | IF(lwp) WRITE(numout,cform_err) |
---|
412 | IF(lwp) WRITE(numout,*)' zgr_bat : unable to read the file', clname |
---|
413 | nstop = nstop + 1 |
---|
414 | ELSE |
---|
415 | zdta(:,:) = 0.e0 |
---|
416 | IF(lwp) WRITE(numout,*)' zgr_bat : bathy_meter not found, but not used, bathy array set to zero' |
---|
417 | ENDIF |
---|
418 | ENDIF |
---|
419 | ! ! =============== ! |
---|
420 | ELSE ! error ! |
---|
421 | ! ! =============== ! |
---|
422 | IF(lwp) WRITE(numout,cform_err) |
---|
423 | IF(lwp) WRITE(numout,*) ' parameter , ntopo = ', ntopo |
---|
424 | nstop = nstop + 1 |
---|
425 | ENDIF |
---|
426 | |
---|
427 | |
---|
428 | ! ======================================= |
---|
429 | ! local domain level and meter bathymetry (mbathy,bathy) |
---|
430 | ! ======================================= |
---|
431 | |
---|
432 | mbathy(:,:) = 0 ! set to zero extra halo points |
---|
433 | bathy (:,:) = 0.e0 ! (require for mpp case) |
---|
434 | |
---|
435 | DO jj = 1, nlcj ! interior values |
---|
436 | DO ji = 1, nlci |
---|
437 | mbathy(ji,jj) = idta( mig(ji), mjg(jj) ) |
---|
438 | bathy (ji,jj) = zdta( mig(ji), mjg(jj) ) |
---|
439 | END DO |
---|
440 | END DO |
---|
441 | |
---|
442 | ! ======================= |
---|
443 | ! NO closed seas or lakes |
---|
444 | ! ======================= |
---|
445 | |
---|
446 | IF( .NOT. lclosea ) THEN |
---|
447 | DO jl = 1, jpncs |
---|
448 | DO jj = ncsj1(jl), ncsj2(jl) |
---|
449 | DO ji = ncsi1(jl), ncsi2(jl) |
---|
450 | mbathy(ji,jj) = 0 ! suppress closed seas |
---|
451 | bathy (ji,jj) = 0.e0 ! and lakes |
---|
452 | END DO |
---|
453 | END DO |
---|
454 | END DO |
---|
455 | ENDIF |
---|
456 | |
---|
457 | ! =========== |
---|
458 | ! Zoom domain |
---|
459 | ! =========== |
---|
460 | |
---|
461 | IF( lzoom ) CALL zgr_bat_zoom |
---|
462 | |
---|
463 | ! ================ |
---|
464 | ! Bathymetry check |
---|
465 | ! ================ |
---|
466 | |
---|
467 | CALL zgr_bat_ctl |
---|
468 | |
---|
469 | END SUBROUTINE zgr_bat |
---|
470 | |
---|
471 | |
---|
472 | SUBROUTINE zgr_bat_zoom |
---|
473 | !!---------------------------------------------------------------------- |
---|
474 | !! *** ROUTINE zgr_bat_zoom *** |
---|
475 | !! |
---|
476 | !! ** Purpose : - Close zoom domain boundary if necessary |
---|
477 | !! - Suppress Med Sea from ORCA R2 and R05 arctic zoom |
---|
478 | !! |
---|
479 | !! ** Method : |
---|
480 | !! |
---|
481 | !! ** Action : - update mbathy: level bathymetry (in level index) |
---|
482 | !! |
---|
483 | !! History : |
---|
484 | !! 9.0 ! 03-08 (G. Madec) Original code |
---|
485 | !!---------------------------------------------------------------------- |
---|
486 | !! * Local variables |
---|
487 | INTEGER :: ii0, ii1, ij0, ij1 ! temporary integers |
---|
488 | !!---------------------------------------------------------------------- |
---|
489 | |
---|
490 | IF(lwp) WRITE(numout,*) |
---|
491 | IF(lwp) WRITE(numout,*) ' zgr_bat_zoom : modify the level bathymetry for zoom domain' |
---|
492 | IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~' |
---|
493 | |
---|
494 | ! Zoom domain |
---|
495 | ! =========== |
---|
496 | |
---|
497 | ! Forced closed boundary if required |
---|
498 | IF( lzoom_w ) mbathy( mi0(jpizoom):mi1(jpizoom) , : ) = 0 |
---|
499 | IF( lzoom_s ) mbathy( : , mj0(jpjzoom):mj1(jpjzoom) ) = 0 |
---|
500 | IF( lzoom_e ) mbathy( mi0(jpiglo+jpizoom-1):mi1(jpiglo+jpizoom-1) , : ) = 0 |
---|
501 | IF( lzoom_n ) mbathy( : , mj0(jpjglo+jpjzoom-1):mj1(jpjglo+jpjzoom-1) ) = 0 |
---|
502 | |
---|
503 | ! Configuration specific domain modifications |
---|
504 | ! (here, ORCA arctic configuration: suppress Med Sea) |
---|
505 | IF( cp_cfg == "orca" .AND. lzoom_arct ) THEN |
---|
506 | SELECT CASE ( jp_cfg ) |
---|
507 | ! ! ======================= |
---|
508 | CASE ( 2 ) ! ORCA_R2 configuration |
---|
509 | ! ! ======================= |
---|
510 | IF(lwp) WRITE(numout,*) ' ORCA R2 arctic zoom: suppress the Med Sea' |
---|
511 | |
---|
512 | ii0 = 141 ; ii1 = 162 ! Sea box i,j indices |
---|
513 | ij0 = 98 ; ij1 = 110 |
---|
514 | ! ! ======================= |
---|
515 | CASE ( 05 ) ! ORCA_R05 configuration |
---|
516 | ! ! ======================= |
---|
517 | IF(lwp) WRITE(numout,*) ' ORCA R05 arctic zoom: suppress the Med Sea' |
---|
518 | |
---|
519 | ii0 = 563 ; ii1 = 642 ! zero over the Med Sea boxe |
---|
520 | ij0 = 314 ; ij1 = 370 |
---|
521 | |
---|
522 | END SELECT |
---|
523 | ! |
---|
524 | mbathy( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0 ! zero over the Med Sea boxe |
---|
525 | ! |
---|
526 | ENDIF |
---|
527 | |
---|
528 | END SUBROUTINE zgr_bat_zoom |
---|
529 | |
---|
530 | |
---|
531 | SUBROUTINE zgr_bat_ctl |
---|
532 | !!---------------------------------------------------------------------- |
---|
533 | !! *** ROUTINE zgr_bat_ctl *** |
---|
534 | !! |
---|
535 | !! ** Purpose : check the bathymetry in levels |
---|
536 | !! |
---|
537 | !! ** Method : The array mbathy is checked to verified its consistency |
---|
538 | !! with the model options. in particular: |
---|
539 | !! mbathy must have at least 1 land grid-points (mbathy<=0) |
---|
540 | !! along closed boundary. |
---|
541 | !! mbathy must be cyclic IF jperio=1. |
---|
542 | !! mbathy must be lower or equal to jpk-1. |
---|
543 | !! isolated ocean grid points are suppressed from mbathy |
---|
544 | !! since they are only connected to remaining |
---|
545 | !! ocean through vertical diffusion. |
---|
546 | !! C A U T I O N : mbathy will be modified during the initializa- |
---|
547 | !! tion phase to become the number of non-zero w-levels of a water |
---|
548 | !! column, with a minimum value of 1. |
---|
549 | !! |
---|
550 | !! ** Action : - update mbathy: level bathymetry (in level index) |
---|
551 | !! - update bathy : meter bathymetry (in meters) |
---|
552 | !! |
---|
553 | !! History : |
---|
554 | !! 9.0 ! 03-08 (G. Madec) Original code |
---|
555 | !!---------------------------------------------------------------------- |
---|
556 | !! * Local declarations |
---|
557 | INTEGER :: ji, jj, jl ! dummy loop indices |
---|
558 | INTEGER :: & |
---|
559 | icompt, ibtest, ikmax ! temporary integers |
---|
560 | REAL(wp), DIMENSION(jpi,jpj) :: & |
---|
561 | zbathy ! temporary workspace |
---|
562 | !!---------------------------------------------------------------------- |
---|
563 | |
---|
564 | IF(lwp) WRITE(numout,*) |
---|
565 | IF(lwp) WRITE(numout,*) ' zgr_bat_ctl : check the bathymetry' |
---|
566 | IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~' |
---|
567 | |
---|
568 | ! ================ |
---|
569 | ! Bathymetry check |
---|
570 | ! ================ |
---|
571 | |
---|
572 | ! Suppress isolated ocean grid points |
---|
573 | |
---|
574 | IF(lwp) WRITE(numout,*) |
---|
575 | IF(lwp) WRITE(numout,*)' suppress isolated ocean grid points' |
---|
576 | IF(lwp) WRITE(numout,*)' -----------------------------------' |
---|
577 | |
---|
578 | icompt = 0 |
---|
579 | |
---|
580 | DO jl = 1, 2 |
---|
581 | |
---|
582 | IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN |
---|
583 | mbathy( 1 ,:) = mbathy(jpim1,:) |
---|
584 | mbathy(jpi,:) = mbathy( 2 ,:) |
---|
585 | ENDIF |
---|
586 | DO jj = 2, jpjm1 |
---|
587 | DO ji = 2, jpim1 |
---|
588 | ibtest = MAX( mbathy(ji-1,jj), mbathy(ji+1,jj), & |
---|
589 | mbathy(ji,jj-1),mbathy(ji,jj+1) ) |
---|
590 | IF( ibtest < mbathy(ji,jj) ) THEN |
---|
591 | IF(lwp) WRITE(numout,*) ' the number of ocean level at ', & |
---|
592 | 'grid-point (i,j) = ',ji,jj,' is changed from ', & |
---|
593 | mbathy(ji,jj),' to ', ibtest |
---|
594 | mbathy(ji,jj) = ibtest |
---|
595 | icompt = icompt + 1 |
---|
596 | ENDIF |
---|
597 | END DO |
---|
598 | END DO |
---|
599 | |
---|
600 | END DO |
---|
601 | IF( icompt == 0 ) THEN |
---|
602 | IF(lwp) WRITE(numout,*)' no isolated ocean grid points' |
---|
603 | ELSE |
---|
604 | IF(lwp) WRITE(numout,*)' ',icompt,' ocean grid points suppressed' |
---|
605 | ENDIF |
---|
606 | IF( lk_mpp ) THEN |
---|
607 | zbathy(:,:) = FLOAT( mbathy(:,:) ) |
---|
608 | CALL lbc_lnk( zbathy, 'T', 1. ) |
---|
609 | mbathy(:,:) = INT( zbathy(:,:) ) |
---|
610 | ENDIF |
---|
611 | |
---|
612 | ! 3.2 East-west cyclic boundary conditions |
---|
613 | |
---|
614 | IF( nperio == 0 ) THEN |
---|
615 | IF(lwp) WRITE(numout,*) ' mbathy set to 0 along east and west', & |
---|
616 | ' boundary: nperio = ', nperio |
---|
617 | IF( lk_mpp ) THEN |
---|
618 | IF( nbondi == -1 .OR. nbondi == 2 ) THEN |
---|
619 | IF( jperio /= 1 ) mbathy(1,:) = 0 |
---|
620 | ENDIF |
---|
621 | IF( nbondi == 1 .OR. nbondi == 2 ) THEN |
---|
622 | IF( jperio /= 1 ) mbathy(nlci,:) = 0 |
---|
623 | ENDIF |
---|
624 | ELSE |
---|
625 | mbathy( 1 ,:) = 0 |
---|
626 | mbathy(jpi,:) = 0 |
---|
627 | ENDIF |
---|
628 | ELSEIF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN |
---|
629 | IF(lwp) WRITE(numout,*)' east-west cyclic boundary conditions', & |
---|
630 | ' on mbathy: nperio = ', nperio |
---|
631 | mbathy( 1 ,:) = mbathy(jpim1,:) |
---|
632 | mbathy(jpi,:) = mbathy( 2 ,:) |
---|
633 | ELSEIF( nperio == 2 ) THEN |
---|
634 | IF(lwp) WRITE(numout,*) ' equatorial boundary conditions', & |
---|
635 | ' on mbathy: nperio = ', nperio |
---|
636 | ELSE |
---|
637 | IF(lwp) WRITE(numout,*) ' e r r o r' |
---|
638 | IF(lwp) WRITE(numout,*) ' parameter , nperio = ', nperio |
---|
639 | ! STOP 'dom_mba' |
---|
640 | ENDIF |
---|
641 | |
---|
642 | ! Set to zero mbathy over islands if necessary (lk_isl=F) |
---|
643 | IF( .NOT. lk_isl ) THEN ! No island |
---|
644 | IF(lwp) WRITE(numout,*) |
---|
645 | IF(lwp) WRITE(numout,*) ' mbathy set to 0 over islands' |
---|
646 | IF(lwp) WRITE(numout,*) ' ----------------------------' |
---|
647 | |
---|
648 | mbathy(:,:) = MAX( 0, mbathy(:,:) ) |
---|
649 | |
---|
650 | ! Boundary condition on mbathy |
---|
651 | IF( .NOT.lk_mpp ) THEN |
---|
652 | !!bug ??? y reflechir! |
---|
653 | ! ... mono- or macro-tasking: T-point, >0, 2D array, no slab |
---|
654 | zbathy(:,:) = FLOAT( mbathy(:,:) ) |
---|
655 | CALL lbc_lnk( zbathy, 'T', 1. ) |
---|
656 | mbathy(:,:) = INT( zbathy(:,:) ) |
---|
657 | ENDIF |
---|
658 | |
---|
659 | ENDIF |
---|
660 | |
---|
661 | ! Number of ocean level inferior or equal to jpkm1 |
---|
662 | |
---|
663 | ikmax = 0 |
---|
664 | DO jj = 1, jpj |
---|
665 | DO ji = 1, jpi |
---|
666 | ikmax = MAX( ikmax, mbathy(ji,jj) ) |
---|
667 | END DO |
---|
668 | END DO |
---|
669 | !!! test a faire: ikmax = MAX( mbathy(:,:) ) ??? |
---|
670 | |
---|
671 | IF( ikmax > jpkm1 ) THEN |
---|
672 | IF(lwp) WRITE(numout,*) ' maximum number of ocean level = ', ikmax,' > jpk-1' |
---|
673 | IF(lwp) WRITE(numout,*) ' change jpk to ',ikmax+1,' to use the exact ead bathymetry' |
---|
674 | ELSE IF( ikmax < jpkm1 ) THEN |
---|
675 | IF(lwp) WRITE(numout,*) ' maximum number of ocean level = ', ikmax,' < jpk-1' |
---|
676 | IF(lwp) WRITE(numout,*) ' you can decrease jpk to ', ikmax+1 |
---|
677 | ENDIF |
---|
678 | |
---|
679 | IF( lwp .AND. nprint == 1 ) THEN |
---|
680 | WRITE(numout,*) |
---|
681 | WRITE(numout,*) ' bathymetric field ' |
---|
682 | WRITE(numout,*) ' ------------------' |
---|
683 | WRITE(numout,*) ' number of non-zero T-levels ' |
---|
684 | CALL prihin( mbathy, jpi, jpj, 1, jpi, & |
---|
685 | 1 , 1 , jpj, 1, 3 , & |
---|
686 | numout ) |
---|
687 | WRITE(numout,*) |
---|
688 | ENDIF |
---|
689 | |
---|
690 | END SUBROUTINE zgr_bat_ctl |
---|
691 | |
---|
692 | |
---|
693 | # include "domzgr_zps.h90" |
---|
694 | |
---|
695 | |
---|
696 | # include "domzgr_s.h90" |
---|
697 | |
---|
698 | |
---|
699 | !!====================================================================== |
---|
700 | END MODULE domzgr |
---|