ROMS
Loading...
Searching...
No Matches
ana_grid.h
Go to the documentation of this file.
1!
2 SUBROUTINE ana_grid (ng, tile, model)
3!
4!! git $Id$
5!!======================================================================
6!! Copyright (c) 2002-2025 The ROMS Group !
7!! Licensed under a MIT/X style license !
8!! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine sets model grid using an analytical expressions. !
12! !
13! On Output: stored in common blocks: !
14! !
15! "grid" (file grid.h) !
16! "scalars" (file scalar.h) !
17! !
18! el Length (m) of domain box in the ETA-direction. !
19! f Coriolis parameter (1/seconds) at RHO-points. !
20! h Bathymetry (meters; positive) at RHO-points. !
21! hmin Minimum depth of bathymetry (m). !
22! hmax Maximum depth of bathymetry (m). !
23! pm Coordinate transformation metric "m" (1/meters) !
24! associated with the differential distances in XI !
25! at RHO-points. !
26! pn Coordinate transformation metric "n" (1/meters) !
27! associated with the differential distances in ETA. !
28! at RHO-points. !
29! xl Length (m) of domain box in the XI-direction. !
30! xp XI-coordinates (m) at PSI-points. !
31! xr XI-coordinates (m) at RHO-points. !
32! yp ETA-coordinates (m) at PSI-points. !
33! yr ETA-coordinates (m) at RHO-points. !
34! !
35!=======================================================================
36!
37 USE mod_param
38 USE mod_grid
39 USE mod_ncparam
40!
41! Imported variable declarations.
42!
43 integer, intent(in) :: ng, tile, model
44!
45! Local variable declarations.
46!
47 character (len=*), parameter :: MyFile = &
48 & __FILE__
49!
50#include "tile.h"
51!
52 CALL ana_grid_tile (ng, tile, model, &
53 & lbi, ubi, lbj, ubj, &
54 & imins, imaxs, jmins, jmaxs, &
55 & grid(ng) % angler, &
56#if defined CURVGRID && defined UV_ADV
57 & grid(ng) % dmde, &
58 & grid(ng) % dndx, &
59#endif
60#ifdef ICESHELF
61 & grid(ng) % zice, &
62#endif
63#ifdef SPHERICAL
64 & grid(ng) % lonp, &
65 & grid(ng) % lonr, &
66 & grid(ng) % lonu, &
67 & grid(ng) % lonv, &
68 & grid(ng) % latp, &
69 & grid(ng) % latr, &
70 & grid(ng) % latu, &
71 & grid(ng) % latv, &
72#else
73 & grid(ng) % xp, &
74 & grid(ng) % xr, &
75 & grid(ng) % xu, &
76 & grid(ng) % xv, &
77 & grid(ng) % yp, &
78 & grid(ng) % yr, &
79 & grid(ng) % yu, &
80 & grid(ng) % yv, &
81#endif
82 & grid(ng) % pn, &
83 & grid(ng) % pm, &
84 & grid(ng) % f, &
85 & grid(ng) % h)
86!
87! Set analytical header file name used.
88!
89#ifdef DISTRIBUTE
90 IF (lanafile) THEN
91#else
92 IF (lanafile.and.(tile.eq.0)) THEN
93#endif
94 ananame( 7)=myfile
95 END IF
96!
97 RETURN
98 END SUBROUTINE ana_grid
99!
100!***********************************************************************
101 SUBROUTINE ana_grid_tile (ng, tile, model, &
102 & LBi, UBi, LBj, UBj, &
103 & IminS, ImaxS, JminS, JmaxS, &
104 & angler, &
105#if defined CURVGRID && defined UV_ADV
106 & dmde, dndx, &
107#endif
108#ifdef ICESHELF
109 & zice, &
110#endif
111#ifdef SPHERICAL
112 & lonp, lonr, lonu, lonv, &
113 & latp, latr, latu, latv, &
114#else
115 & xp, xr, xu, xv, &
116 & yp, yr, yu, yv, &
117#endif
118 & pn, pm, f, h)
119!***********************************************************************
120!
121 USE mod_param
122 USE mod_parallel
123 USE mod_ncparam
124 USE mod_iounits
125 USE mod_scalars
126!
128#ifdef DISTRIBUTE
130#endif
131 USE stats_mod, ONLY : stats_2dfld
132!
133! Imported variable declarations.
134!
135 integer, intent(in) :: ng, tile, model
136 integer, intent(in) :: LBi, UBi, LBj, UBj
137 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
138!
139#ifdef ASSUMED_SHAPE
140 real(r8), intent(out) :: angler(LBi:,LBj:)
141# if defined CURVGRID && defined UV_ADV
142 real(r8), intent(out) :: dmde(LBi:,LBj:)
143 real(r8), intent(out) :: dndx(LBi:,LBj:)
144# endif
145# ifdef ICESHELF
146 real(r8), intent(out) :: zice(LBi:,LBj:)
147# endif
148# ifdef SPHERICAL
149 real(r8), intent(out) :: lonp(LBi:,LBj:)
150 real(r8), intent(out) :: lonr(LBi:,LBj:)
151 real(r8), intent(out) :: lonu(LBi:,LBj:)
152 real(r8), intent(out) :: lonv(LBi:,LBj:)
153 real(r8), intent(out) :: latp(LBi:,LBj:)
154 real(r8), intent(out) :: latr(LBi:,LBj:)
155 real(r8), intent(out) :: latu(LBi:,LBj:)
156 real(r8), intent(out) :: latv(LBi:,LBj:)
157# else
158 real(r8), intent(out) :: xp(LBi:,LBj:)
159 real(r8), intent(out) :: xr(LBi:,LBj:)
160 real(r8), intent(out) :: xu(LBi:,LBj:)
161 real(r8), intent(out) :: xv(LBi:,LBj:)
162 real(r8), intent(out) :: yp(LBi:,LBj:)
163 real(r8), intent(out) :: yr(LBi:,LBj:)
164 real(r8), intent(out) :: yu(LBi:,LBj:)
165 real(r8), intent(out) :: yv(LBi:,LBj:)
166# endif
167 real(r8), intent(out) :: pn(LBi:,LBj:)
168 real(r8), intent(out) :: pm(LBi:,LBj:)
169 real(r8), intent(out) :: f(LBi:,LBj:)
170 real(r8), intent(out) :: h(LBi:,LBj:)
171#else
172 real(r8), intent(out) :: angler(LBi:UBi,LBj:UBj)
173# if defined CURVGRID && defined UV_ADV
174 real(r8), intent(out) :: dmde(LBi:UBi,LBj:UBj)
175 real(r8), intent(out) :: dndx(LBi:UBi,LBj:UBj)
176# endif
177# ifdef ICESHELF
178 real(r8), intent(out) :: zice(LBi:UBi,LBj:UBj)
179# endif
180# ifdef SPHERICAL
181 real(r8), intent(out) :: lonp(LBi:UBi,LBj:UBj)
182 real(r8), intent(out) :: lonr(LBi:UBi,LBj:UBj)
183 real(r8), intent(out) :: lonu(LBi:UBi,LBj:UBj)
184 real(r8), intent(out) :: lonv(LBi:UBi,LBj:UBj)
185 real(r8), intent(out) :: latp(LBi:UBi,LBj:UBj)
186 real(r8), intent(out) :: latr(LBi:UBi,LBj:UBj)
187 real(r8), intent(out) :: latu(LBi:UBi,LBj:UBj)
188 real(r8), intent(out) :: latv(LBi:UBi,LBj:UBj)
189# else
190 real(r8), intent(out) :: xp(LBi:UBi,LBj:UBj)
191 real(r8), intent(out) :: xr(LBi:UBi,LBj:UBj)
192 real(r8), intent(out) :: xu(LBi:UBi,LBj:UBj)
193 real(r8), intent(out) :: xv(LBi:UBi,LBj:UBj)
194 real(r8), intent(out) :: yp(LBi:UBi,LBj:UBj)
195 real(r8), intent(out) :: yr(LBi:UBi,LBj:UBj)
196 real(r8), intent(out) :: yu(LBi:UBi,LBj:UBj)
197 real(r8), intent(out) :: yv(LBi:UBi,LBj:UBj)
198# endif
199 real(r8), intent(out) :: pn(LBi:UBi,LBj:UBj)
200 real(r8), intent(out) :: pm(LBi:UBi,LBj:UBj)
201 real(r8), intent(out) :: f(LBi:UBi,LBj:UBj)
202 real(r8), intent(out) :: h(LBi:UBi,LBj:UBj)
203#endif
204!
205! Local variable declarations.
206!
207 logical, save :: first = .true.
208!
209 integer :: Imin, Imax, Jmin, Jmax
210 integer :: i, ival, j, k
211!
212 real(r8), parameter :: twopi = 2.0_r8*pi
213
214 real(r8) :: Esize, Xsize, beta, cff, depth, dth
215 real(r8) :: dx, dy, f0, r, theta, val1, val2
216
217#ifdef WEDDELL
218 real(r8) :: hwrk(-1:235), xwrk(-1:235), zwrk
219#endif
220 real(r8) :: wrkX(IminS:ImaxS,JminS:JmaxS)
221 real(r8) :: wrkY(IminS:ImaxS,JminS:JmaxS)
222!
223 TYPE (T_STATS), save :: Stats(16)
224
225#include "set_bounds.h"
226!
227!-----------------------------------------------------------------------
228! Set grid parameters:
229!
230! Xsize Length (m) of domain box in the XI-direction.
231! Esize Length (m) of domain box in the ETA-direction.
232! depth Maximum depth of bathymetry (m).
233! f0 Coriolis parameter, f-plane constant (1/s).
234! beta Coriolis parameter, beta-plane constant (1/s/m).
235!-----------------------------------------------------------------------
236!
237#if defined BASIN
238 xsize=3600.0e+03_r8
239 esize=2800.0e+03_r8
240 depth=5000.0_r8
241 f0=1.0e-04_r8
242 beta=2.0e-11_r8
243#elif defined BENCHMARK
244 xsize=360.0_r8 ! degrees of longitude
245 esize=20.0_r8 ! degrees of latitude
246 depth=4000.0_r8
247 f0=-1.0e-04_r8
248 beta=2.0e-11_r8
249#elif defined BL_TEST
250 xsize=100.0e+03_r8
251 esize=5.0e+03_r8
252 depth=47.5_r8
253 f0=9.25e-04_r8
254 beta=0.0_r8
255#elif defined CHANNEL
256 xsize=600.0e+03_r8
257 esize=360.0e+03_r8
258 depth=500.0_r8
259 f0=1.0e-04_r8
260 beta=0.0_r8
261#elif defined CANYON
262 xsize=128.0e+03_r8
263 esize=96.0e+03_r8
264 depth=4000.0_r8
265 f0=1.0e-04_r8
266 beta=0.0_r8
267#elif defined COUPLING_TEST
268 xsize=6000.0_r8*real(lm(ng),r8)
269 esize=6000.0_r8*real(mm(ng),r8)
270 depth=1500.0_r8
271 f0=5.0e-05_r8
272 beta=0.0_r8
273#elif defined DOUBLE_GYRE
274 xsize=1000.0e+03_r8
275 esize=2000.0e+03_r8
276 depth=500.0_r8
277!! depth=5000.0_r8
278 f0=7.3e-05_r8
279 beta=2.0e-11_r8
280#elif defined ESTUARY_TEST
281 xsize=100000.0_r8
282 esize=300.0_r8
283 depth=10.0_r8
284 f0=0.0_r8
285 beta=0.0_r8
286#elif defined KELVIN
287 xsize=20000.0_r8*real(lm(ng),r8)
288 esize=20000.0_r8*real(mm(ng),r8)
289 depth=100.0_r8
290 f0=1.0e-04_r8
291 beta=0.0_r8
292#elif defined FLT_TEST
293 xsize=1.0e+03_r8*real(lm(ng),r8)
294 esize=1.0e+03_r8*real(mm(ng),r8)
295 depth=10.0_r8
296 f0=0.0_r8
297 beta=0.0_r8
298#elif defined GRAV_ADJ
299 xsize=64.0e+03_r8
300!! Esize=2.0E+03_r8
301 esize=mm(ng)*xsize/lm(ng)
302 depth=20.0_r8
303 f0=0.0_r8
304 beta=0.0_r8
305#elif defined LAB_CANYON
306 xsize=0.55_r8 ! width of annulus
307 esize=2.0_r8*pi ! azimuthal length (radians)
308 f0=4.0_r8*pi/25.0_r8
309 beta=0.0_r8
310#elif defined LAKE_SIGNELL
311 xsize=50.0e3_r8
312 esize=10.0e3_r8
313 depth=18.0_r8
314 f0=0.0e-04_r8
315 beta=0.0_r8
316#elif defined LMD_TEST
317 xsize=100.0e+03_r8
318 esize=100.0e+03_r8
319 depth=50.0_r8
320 f0=1.09e-04_r8
321 beta=0.0_r8
322# elif defined MIXED_LAYER
323 xsize=500.0_r8
324 esize=400.0_r8
325 depth=50.0_r8
326 f0=0.0_r8
327 beta=0.0_r8
328#elif defined OVERFLOW
329 xsize=4.0e+03_r8
330 esize=200.0e+03_r8
331 depth=4000.0_r8
332 f0=0.0_r8
333 beta=0.0_r8
334#elif defined RIVERPLUME1
335 xsize=58.5e+03_r8
336 esize=201.0e+03_r8
337 depth=150.0_r8
338 f0=1.0e-04_r8
339 beta=0.0_r8
340#elif defined RIVERPLUME2
341 xsize=100.0e+03_r8
342 esize=210.0e+03_r8
343 depth=190.0_r8
344 f0=1.0e-04_r8
345 beta=0.0_r8
346#elif defined SEAMOUNT
347 xsize=320.0e+03_r8
348 esize=320.0e+03_r8
349 depth=5000.0_r8
350!! f0=1.0E-04_r8
351 f0=0.0_r8
352 beta=0.0_r8
353#elif defined SOLITON
354!! Xsize=0.5_r8*REAL(Lm(ng),r8)
355!! Esize=0.5_r8*REAL(Mm(ng),r8)
356 xsize=48.0_r8
357 esize=16.0_r8
358 depth=1.0_r8
359 f0=0.0_r8
360 beta=1.0_r8
361 g=1.0_r8
362#elif defined SED_TEST1
363 xsize=300.0_r8
364 esize=36.0_r8
365 depth=10.0_r8
366 f0=0.0_r8
367 beta=0.0_r8
368#elif defined SED_TOY
369 xsize=40.0_r8
370 esize=30.0_r8
371 depth=0.5_r8
372 f0=0.0_r8
373 beta=0.0_r8
374# elif defined SHOREFACE
375 xsize=1180.0_r8
376 esize=140.0_r8
377 depth=15.0_r8
378 f0=0.0e-04_r8
379 beta=0.0_r8
380#elif defined TEST_CHAN
381 xsize=10000.0_r8
382 esize=1000.0_r8
383 depth=10.0_r8
384 f0=0.0_r8
385 beta=0.0_r8
386#elif defined UPWELLING
387 xsize=1000.0_r8*real(lm(ng),r8)
388 esize=1000.0_r8*real(mm(ng),r8)
389 depth=150.0_r8
390 f0=-8.26e-05_r8
391 beta=0.0_r8
392#elif defined WEDDELL
393 xsize=4000.0_r8*real(lm(ng),r8)
394 esize=4000.0_r8*real(mm(ng),r8)
395 depth=4500.0_r8
396 f0=0.0_r8
397 beta=0.0_r8
398#elif defined WINDBASIN
399 xsize=2000.0_r8*real(lm(ng),r8)
400 esize=1000.0_r8*real(mm(ng),r8)
401 depth=50.0_r8
402 f0=1.0e-04_r8
403 beta=0.0_r8
404#else
405 ana_grid.h: no values provided for xsize, esize, depth, f0, beta.
406#endif
407!
408! Load grid parameters to global storage.
409!
410 IF (domain(ng)%NorthEast_Test(tile)) THEN
411 xl(ng)=xsize
412 el(ng)=esize
413 END IF
414!
415!-----------------------------------------------------------------------
416! Initialize field statistics structure.
417!-----------------------------------------------------------------------
418!
419 IF (first) THEN
420 first=.false.
421 DO i=1,SIZE(stats,1)
422 stats(i) % checksum=0_i8b
423 stats(i) % count=0
424 stats(i) % min=large
425 stats(i) % max=-large
426 stats(i) % avg=0.0_r8
427 stats(i) % rms=0.0_r8
428 END DO
429 END IF
430 IF (domain(ng)%NorthEast_Corner(tile)) WRITE (stdout,'(1x)')
431!
432!-----------------------------------------------------------------------
433! Compute the (XI,ETA) coordinates at PSI- and RHO-points.
434! Set grid spacing (m).
435!-----------------------------------------------------------------------
436!
437! Determine I- and J-ranges for computing grid data. These ranges
438! are special in periodic boundary conditons since periodicity cannot
439! be imposed in the grid coordinates.
440!
441 IF (domain(ng)%Western_Edge(tile)) THEN
442 imin=istr-1
443 ELSE
444 imin=istr
445 END IF
446 IF (domain(ng)%Eastern_Edge(tile)) THEN
447 imax=iend+1
448 ELSE
449 imax=iend
450 END IF
451 IF (domain(ng)%Southern_Edge(tile)) THEN
452 jmin=jstr-1
453 ELSE
454 jmin=jstr
455 END IF
456 IF (domain(ng)%Northern_Edge(tile)) THEN
457 jmax=jend+1
458 ELSE
459 jmax=jend
460 END IF
461
462#if defined BENCHMARK
463!
464! Spherical coordinates set-up.
465!
466 dx=xsize/real(lm(ng),r8)
467 dy=esize/real(mm(ng),r8)
468 spherical=.true.
469 DO j=jmin,jmax
470 val1=-70.0_r8+dy*(real(j,r8)-0.5_r8)
471 val2=-70.0_r8+dy*real(j,r8)
472 DO i=imin,imax
473 lonr(i,j)=dx*(real(i,r8)-0.5_r8)
474 latr(i,j)=val1
475 lonu(i,j)=dx*real(i,r8)
476 lonp(i,j)=lonu(i,j)
477 latu(i,j)=latr(i,j)
478 lonv(i,j)=lonr(i,j)
479 latv(i,j)=val2
480 latp(i,j)=latv(i,j)
481 END DO
482 END DO
483#elif defined LAB_CANYON
484!
485! Polar coordinates set-up.
486!
487 dx=xsize/real(lm(ng),r8)
488 dy=esize/real(mm(ng),r8)
489!! dth=twopi/REAL(Mm(ng),r8) ! equal azimultal spacing
490 dth=0.01_r8 ! azimultal spacing
491 cff=(4.0_r8*pi/(dth*real(mm(ng),r8)))-1.0_r8 ! F
492 DO j=jmin,jmax
493 DO i=imin,imax
494 r=0.35_r8+dx*real(i-1,r8)
495 theta=-pi+ &
496 & 0.5_r8*dth*((cff+1.0_r8)*real(j-1,r8)+ &
497 & (cff-1.0_r8)*(real(mm(ng),r8)/twopi)* &
498 & sin(twopi*real(j-1,r8)/real(mm(ng),r8)))
499 xp(i,j)=r*cos(theta)
500 yp(i,j)=r*sin(theta)
501 r=0.35_r8+dx*(real(i-1,r8)+0.5_r8)
502 theta=-pi+ &
503 & 0.5_r8*dth*((cff+1.0_r8)*(real(j-1,r8)+0.5_r8)+ &
504 & (cff-1.0_r8)*(real(mm(ng),r8)/twopi)* &
505 & sin(twopi*(real(j-1,r8)+0.5_r8)/ &
506 & real(mm(ng),r8)))
507 xr(i,j)=r*cos(theta)
508 yr(i,j)=r*sin(theta)
509 xu(i,j)=xp(i,j)
510 yu(i,j)=yr(i,j)
511 xv(i,j)=xr(i,j)
512 yv(i,j)=yp(i,j)
513 END DO
514 END DO
515#else
516 dx=xsize/real(lm(ng),r8)
517 dy=esize/real(mm(ng),r8)
518 DO j=jmin,jmax
519 DO i=imin,imax
520# ifdef BL_TEST
521 dx=0.5_r8*(4000.0_r8/real(lm(ng)+1,r8))*real(i,r8)+675.0_r8
522# endif
523 xp(i,j)=dx*real(i-1,r8)
524 xr(i,j)=dx*(real(i-1,r8)+0.5_r8)
525 xu(i,j)=xp(i,j)
526 xv(i,j)=xr(i,j)
527 yp(i,j)=dy*real(j-1,r8)
528 yr(i,j)=dy*(real(j-1,r8)+0.5_r8)
529 yu(i,j)=yr(i,j)
530 yv(i,j)=yp(i,j)
531 END DO
532 END DO
533#endif
534!
535! Report statistics.
536!
537#ifdef SPHERICAL
538 CALL stats_2dfld (ng, tile, inlm, p2dvar, stats(1), 0, &
539 & lbi, ubi, lbj, ubj, lonp)
540 IF (domain(ng)%NorthEast_Corner(tile)) THEN
541 WRITE (stdout,10) 'longitude of PSI-points: lon_psi', &
542 & ng, stats(1)%min, stats(1)%max
543 END IF
544 CALL stats_2dfld (ng, tile, inlm, p2dvar, stats(2), 0, &
545 & lbi, ubi, lbj, ubj, latp)
546 IF (domain(ng)%NorthEast_Corner(tile)) THEN
547 WRITE (stdout,10) 'latitude of PSI-points: lat_psi', &
548 & ng, stats(2)%min, stats(2)%max
549 END IF
550
551 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(3), 0, &
552 & lbi, ubi, lbj, ubj, lonr)
553 IF (domain(ng)%NorthEast_Corner(tile)) THEN
554 WRITE (stdout,10) 'longitude of RHO-points: lon_rho', &
555 & ng, stats(3)%min, stats(3)%max
556 END IF
557 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(4), 0, &
558 & lbi, ubi, lbj, ubj, latr)
559 IF (domain(ng)%NorthEast_Corner(tile)) THEN
560 WRITE (stdout,10) 'latitude of RHO-points: lat_rho', &
561 & ng, stats(4)%min, stats(4)%max
562 END IF
563
564 CALL stats_2dfld (ng, tile, inlm, u2dvar, stats(5), 0, &
565 & lbi, ubi, lbj, ubj, lonu)
566 IF (domain(ng)%NorthEast_Corner(tile)) THEN
567 WRITE (stdout,10) 'longitude of U-points: lon_u', &
568 & ng, stats(5)%min, stats(5)%max
569 END IF
570 CALL stats_2dfld (ng, tile, inlm, u2dvar, stats(6), 0, &
571 & lbi, ubi, lbj, ubj, latu)
572 IF (domain(ng)%NorthEast_Corner(tile)) THEN
573 WRITE (stdout,10) 'latitude of U-points: lat_u', &
574 & ng, stats(6)%min, stats(6)%max
575 END IF
576
577 CALL stats_2dfld (ng, tile, inlm, v2dvar, stats(7), 0, &
578 & lbi, ubi, lbj, ubj, lonv)
579 IF (domain(ng)%NorthEast_Corner(tile)) THEN
580 WRITE (stdout,10) 'longitude of V-points: lon_v', &
581 & ng, stats(7)%min, stats(7)%max
582 END IF
583 CALL stats_2dfld (ng, tile, inlm, v2dvar, stats(8), 0, &
584 & lbi, ubi, lbj, ubj, latv)
585 IF (domain(ng)%NorthEast_Corner(tile)) THEN
586 WRITE (stdout,10) 'latitude of V-points: lat_v', &
587 & ng, stats(8)%min, stats(8)%max
588 END IF
589#else
590 CALL stats_2dfld (ng, tile, inlm, p2dvar, stats(1), 0, &
591 & lbi, ubi, lbj, ubj, xp)
592 IF (domain(ng)%NorthEast_Corner(tile)) THEN
593 WRITE (stdout,10) 'x-location of PSI-points: x_psi', &
594 & ng, stats(1)%min, stats(1)%max
595 END IF
596 CALL stats_2dfld (ng, tile, inlm, p2dvar, stats(2), 0, &
597 & lbi, ubi, lbj, ubj, yp)
598 IF (domain(ng)%NorthEast_Corner(tile)) THEN
599 WRITE (stdout,10) 'y-location of PSI-points: y_psi', &
600 & ng, stats(2)%min, stats(2)%max
601 END IF
602
603 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(3), 0, &
604 & lbi, ubi, lbj, ubj, xr)
605 IF (domain(ng)%NorthEast_Corner(tile)) THEN
606 WRITE (stdout,10) 'x-location of RHO-points: x_rho', &
607 & ng, stats(3)%min, stats(3)%max
608 END IF
609 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(4), 0, &
610 & lbi, ubi, lbj, ubj, yr)
611 IF (domain(ng)%NorthEast_Corner(tile)) THEN
612 WRITE (stdout,10) 'y-location of RHO-points: y_rho', &
613 & ng, stats(4)%min, stats(4)%max
614 END IF
615
616 CALL stats_2dfld (ng, tile, inlm, u2dvar, stats(5), 0, &
617 & lbi, ubi, lbj, ubj, xu)
618 IF (domain(ng)%NorthEast_Corner(tile)) THEN
619 WRITE (stdout,10) 'x-location of U-points: x_u', &
620 & ng, stats(5)%min, stats(5)%max
621 END IF
622 CALL stats_2dfld (ng, tile, inlm, u2dvar, stats(6), 0, &
623 & lbi, ubi, lbj, ubj, yu)
624 IF (domain(ng)%NorthEast_Corner(tile)) THEN
625 WRITE (stdout,10) 'y-location of U-points: y_u', &
626 & ng, stats(6)%min, stats(6)%max
627 END IF
628
629 CALL stats_2dfld (ng, tile, inlm, v2dvar, stats(7), 0, &
630 & lbi, ubi, lbj, ubj, xv)
631 IF (domain(ng)%NorthEast_Corner(tile)) THEN
632 WRITE (stdout,10) 'x-location of V-points: x_v', &
633 & ng, stats(7)%min, stats(7)%max
634 END IF
635 CALL stats_2dfld (ng, tile, inlm, v2dvar, stats(8), 0, &
636 & lbi, ubi, lbj, ubj, yv)
637 IF (domain(ng)%NorthEast_Corner(tile)) THEN
638 WRITE (stdout,10) 'y-location of V-points: y_v', &
639 & ng, stats(8)%min, stats(8)%max
640 END IF
641#endif
642
643#ifdef DISTRIBUTE
644!
645! Exchange boundary data.
646!
647# ifdef SPHERICAL
648 CALL mp_exchange2d (ng, tile, model, 4, &
649 & lbi, ubi, lbj, ubj, &
650 & nghostpoints, .false., .false., &
651 & lonp, lonr, lonu, lonv)
652 CALL mp_exchange2d (ng, tile, model, 4, &
653 & lbi, ubi, lbj, ubj, &
654 & nghostpoints, .false., .false., &
655 & latp, latr, latu, latv)
656# else
657 CALL mp_exchange2d (ng, tile, model, 4, &
658 & lbi, ubi, lbj, ubj, &
659 & nghostpoints, .false., .false., &
660 & xp, xr, xu, xv)
661 CALL mp_exchange2d (ng, tile, model, 4, &
662 & lbi, ubi, lbj, ubj, &
663 & nghostpoints, .false., .false., &
664 & yp, yr, yu, yv)
665# endif
666#endif
667!
668!-----------------------------------------------------------------------
669! Compute coordinate transformation metrics at RHO-points "pm" and
670! "pn" (1/m) associated with the differential distances in XI and
671! ETA, respectively.
672!-----------------------------------------------------------------------
673!
674#define J_RANGE MIN(JstrT,Jstr-1),MAX(Jend+1,JendT)
675#define I_RANGE MIN(IstrT,Istr-1),MAX(Iend+1,IendT)
676
677#if defined BENCHMARK
678!
679! Spherical coordinates set-up.
680!
681 val1=real(lm(ng),r8)/(2.0_r8*pi*eradius)
682 val2=real(mm(ng),r8)*360.0_r8/(2.0_r8*pi*eradius*esize)
683 DO j=j_range
684 cff=1.0_r8/cos((-70.0_r8+dy*(real(j,r8)-0.5_r8))*deg2rad)
685 DO i=i_range
686 wrkx(i,j)=val1*cff
687 wrky(i,j)=val2
688 END DO
689 END DO
690#elif defined LAB_CANYON
691!
692! Polar coordinates set-up.
693!
694 DO j=j_range
695 DO i=i_range
696 r=0.35_r8+dx*(real(i-1,r8)+0.5_r8)
697 theta=0.5_r8*dth*((cff+1.0_r8)+ &
698 & (cff-1.0_r8)* &
699 & cos(twopi*real(j-1,r8)/real(mm(ng),r8)))
700 wrkx(i,j)=1.0_r8/dx
701 wrky(i,j)=1.0_r8/(r*theta)
702 END DO
703 END DO
704#else
705 DO j=j_range
706 DO i=i_range
707# ifdef BL_TEST
708 dx=0.5_r8*(4000.0_r8/real(lm(ng)+1,r8))*real(i,r8)+675.0_r8
709# endif
710 wrkx(i,j)=1.0_r8/dx
711 wrky(i,j)=1.0_r8/dy
712 END DO
713 END DO
714#endif
715#undef J_RANGE
716#undef I_RANGE
717 DO j=jstrt,jendt
718 DO i=istrt,iendt
719 pm(i,j)=wrkx(i,j)
720 pn(i,j)=wrky(i,j)
721 END DO
722 END DO
723!
724! Report statistics.
725!
726 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(9), 0, &
727 & lbi, ubi, lbj, ubj, pm)
728 IF (domain(ng)%NorthEast_Corner(tile)) THEN
729 WRITE (stdout,10) 'reciprocal XI-grid spacing: pm', &
730 & ng, stats(9)%min, stats(9)%max
731 END IF
732 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(10), 0, &
733 & lbi, ubi, lbj, ubj, pn)
734 IF (domain(ng)%NorthEast_Corner(tile)) THEN
735 WRITE (stdout,10) 'reciprocal ETA-grid spacing: pn', &
736 & ng, stats(10)%min, stats(10)%max
737 END IF
738!
739! Exchange boundary data.
740!
741 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
742 CALL exchange_r2d_tile (ng, tile, &
743 & lbi, ubi, lbj, ubj, &
744 & pm)
745 CALL exchange_r2d_tile (ng, tile, &
746 & lbi, ubi, lbj, ubj, &
747 & pn)
748 END IF
749
750#ifdef DISTRIBUTE
751 CALL mp_exchange2d (ng, tile, model, 2, &
752 & lbi, ubi, lbj, ubj, &
753 & nghostpoints, &
754 & ewperiodic(ng), nsperiodic(ng), &
755 & pm, pn)
756#endif
757
758#if (defined CURVGRID && defined UV_ADV)
759!
760!-----------------------------------------------------------------------
761! Compute d(1/n)/d(xi) and d(1/m)/d(eta) at RHO-points.
762!-----------------------------------------------------------------------
763!
764 DO j=jstr,jend
765 DO i=istr,iend
766 dndx(i,j)=0.5_r8*((1.0_r8/wrky(i+1,j ))- &
767 & (1.0_r8/wrky(i-1,j )))
768 dmde(i,j)=0.5_r8*((1.0_r8/wrkx(i ,j+1))- &
769 & (1.0_r8/wrkx(i ,j-1)))
770 END DO
771 END DO
772!
773! Report statistics.
774!
775 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(11), 0, &
776 & lbi, ubi, lbj, ubj, dmde)
777 IF (domain(ng)%NorthEast_Corner(tile)) THEN
778 WRITE (stdout,10) 'ETA-derivative of inverse metric '// &
779 & 'factor pm: dmde', &
780 & ng, stats(11)%min, stats(11)%max
781 END IF
782 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(12), 0, &
783 & lbi, ubi, lbj, ubj, dndx)
784 IF (domain(ng)%NorthEast_Corner(tile)) THEN
785 WRITE (stdout,10) 'XI-derivative of inverse metric '// &
786 & 'factor pn: dndx', &
787 & ng, stats(12)%min, stats(12)%max
788 END IF
789!
790! Exchange boundary data.
791!
792 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
793 CALL exchange_r2d_tile (ng, tile, &
794 & lbi, ubi, lbj, ubj, &
795 & dndx)
796 CALL exchange_r2d_tile (ng, tile, &
797 & lbi, ubi, lbj, ubj, &
798 & dmde)
799 END IF
800
801# ifdef DISTRIBUTE
802 CALL mp_exchange2d (ng, tile, model, 2, &
803 & lbi, ubi, lbj, ubj, &
804 & nghostpoints, &
805 & ewperiodic(ng), nsperiodic(ng), &
806 & dndx, dmde)
807# endif
808#endif
809!
810!-----------------------------------------------------------------------
811! Angle (radians) between XI-axis and true EAST at RHO-points.
812!-----------------------------------------------------------------------
813!
814#if defined LAB_CANYON
815 DO j=jstrt,jendt
816 DO i=istrt,iendt
817 theta=-pi+ &
818 & 0.5_r8*dth*((cff+1.0_r8)*(real(j-1,r8)+0.5_r8)+ &
819 & (cff-1.0_r8)*(real(mm(ng),r8)/twopi)* &
820 & sin(twopi*(real(j-1,r8)+0.5_r8)/ &
821 & real(mm(ng),r8)))
822 angler(i,j)=theta
823 END DO
824 END DO
825#elif defined WEDDELL
826 val1=90.0_r8*deg2rad
827 DO j=jstrt,jendt
828 DO i=istrt,iendt
829 angler(i,j)=val1
830 END DO
831 END DO
832#else
833 DO j=jstrt,jendt
834 DO i=istrt,iendt
835 angler(i,j)=0.0_r8
836 END DO
837 END DO
838#endif
839!
840! Report Statistics.
841!
842 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(13), 0, &
843 & lbi, ubi, lbj, ubj, angler)
844 IF (domain(ng)%NorthEast_Corner(tile)) THEN
845 WRITE (stdout,10) 'angle between XI-axis and EAST: '// &
846 & 'angler', &
847 & ng, stats(13)%min, stats(13)%max
848 END IF
849!
850! Exchange boundary data.
851!
852 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
853 CALL exchange_r2d_tile (ng, tile, &
854 & lbi, ubi, lbj, ubj, &
855 & angler)
856 END IF
857
858#ifdef DISTRIBUTE
859 CALL mp_exchange2d (ng, tile, model, 1, &
860 & lbi, ubi, lbj, ubj, &
861 & nghostpoints, &
862 & ewperiodic(ng), nsperiodic(ng), &
863 & angler)
864#endif
865!
866!-----------------------------------------------------------------------
867! Compute Coriolis parameter (1/s) at RHO-points.
868!-----------------------------------------------------------------------
869!
870#if defined BENCHMARK
871 val1=2.0_r8*(2.0_r8*pi*366.25_r8/365.25_r8)/86400.0_r8
872 DO j=jstrt,jendt
873 DO i=istrt,iendt
874 f(i,j)=val1*sin(latr(i,j)*deg2rad)
875 END DO
876 END DO
877#elif defined WEDDELL
878 val1=10.4_r8/real(lm(ng),r8)
879 DO j=jstrt,jendt
880 DO i=istrt,iendt
881 f(i,j)=2.0_r8*7.2e-05_r8* &
882 & sin((-79.0_r8+real(i-1,r8)*val1)*deg2rad)
883 END DO
884 END DO
885#else
886 IF (beta.eq.0.0_r8) THEN
887 DO j=jstrt,jendt
888 DO i=istrt,iendt
889 f(i,j)=f0
890 END DO
891 END DO
892 ELSE
893 val1=0.5_r8*esize
894 DO j=jstrt,jendt
895 DO i=istrt,iendt
896 f(i,j)=f0+beta*(yr(i,j)-val1)
897 END DO
898 END DO
899 END IF
900#endif
901!
902! Report Statistics.
903!
904 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(14), 0, &
905 & lbi, ubi, lbj, ubj, f)
906 IF (domain(ng)%NorthEast_Corner(tile)) THEN
907 WRITE (stdout,10) 'Coriolis parameter at RHO-points: f', &
908 & ng, stats(14)%min, stats(14)%max
909 END IF
910!
911! Exchange boundary data.
912!
913 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
914 CALL exchange_r2d_tile (ng, tile, &
915 & lbi, ubi, lbj, ubj, &
916 & f)
917 END IF
918
919#ifdef DISTRIBUTE
920 CALL mp_exchange2d (ng, tile, model, 1, &
921 & lbi, ubi, lbj, ubj, &
922 & nghostpoints, &
923 & ewperiodic(ng), nsperiodic(ng), &
924 & f)
925#endif
926!
927!-----------------------------------------------------------------------
928! Set bathymetry (meters; positive) at RHO-points.
929!-----------------------------------------------------------------------
930!
931#if defined BENCHMARK
932 DO j=jstrt,jendt
933 DO i=istrt,iendt
934 h(i,j)=500.0_r8+1750.0_r8*(1.0+tanh((68.0_r8+latr(i,j))/dy))
935 END DO
936 END DO
937#elif defined BL_TEST
938 DO j=jstrt,jendt
939 DO i=istrt,iendt
940 val1=(xr(i,j)+500.0_r8)/15000.0_r8
941 h(i,j)=14.0_r8+ &
942 & 25.0_r8*(1.0_r8-exp(-pi*xr(i,j)*1.0e-05_r8))- &
943 & 8.0_r8*exp(-val1*val1)
944 END DO
945 END DO
946#elif defined CANYON
947 DO j=jstrt,jendt
948 DO i=istrt,iendt
949 val1=32000.0_r8-16000.0_r8*(sin(pi*xr(i,j)/xsize))**24
950 h(i,j)=20.0_r8+0.5_r8*(depth-20.0_r8)* &
951 & (1.0_r8+tanh((yr(i,j)-val1)/10000.0_r8))
952 END DO
953 END DO
954#elif defined ESTUARY_TEST
955 DO j=jstrt,jendt
956 DO i=istrt,iendt
957 h(i,j)=5.0_r8+(xsize-xr(i,j))/xsize*5.0_r8
958 END DO
959 END DO
960#elif defined LAB_CANYON
961 DO j=jstrt,jendt
962 DO i=istrt,iendt
963 r=0.35_r8+dx*(real(i-1,r8)+0.5_r8)
964 theta=-pi+ &
965 & 0.5_r8*dth*((cff+1.0_r8)*(real(j-1,r8)+0.5_r8)+ &
966 & (cff-1.0_r8)*(real(mm(ng),r8)/twopi)* &
967 & sin(dth*(real(j-1,r8)+0.5_r8)/ &
968 & real(mm(ng),r8)))
969 val1=0.55_r8-0.15_r8*(cos(pi*theta*0.55_r8/0.2_r8)**2) !r_small
970 val2=0.15_r8+0.15_r8*(cos(pi*theta*0.55_r8/0.2_r8)**2) !lambda
971 IF (abs(theta).ge.0.181818181818_r8) THEN
972 IF (r.le.0.55_r8) THEN
973 h(i,j)=0.025_r8 ! shelf
974 ELSE IF (r.ge.0.7_r8) THEN
975 h(i,j)=0.125_r8 ! deep
976 ELSE
977 h(i,j)=0.125_r8-0.1_r8* &
978 & (cos(0.5_r8*pi*(r-0.55_r8)/0.15_r8)**2)
979 END IF
980 ELSE
981 IF (r.le.val1) THEN
982 h(i,j)=0.025_r8 ! shelf
983 ELSE IF (r.ge.0.7_r8) THEN
984 h(i,j)=0.125_r8 ! deep
985 ELSE
986 h(i,j)=0.125_r8-0.1_r8* &
987 & (cos(0.5_r8*pi*(r-val1)/val2)**2)
988 END IF
989 END IF
990 END DO
991 END DO
992#elif defined LAKE_SIGNELL
993 DO j=jstrt,jendt
994 DO i=istrt,iendt
995 h(i,j)=18.0_r8-16.0_r8*real(mm(ng)-j,r8)/real(mm(ng)-1,r8)
996 END DO
997 END DO
998# elif defined MIXED_LAYER
999 DO j=jstrt,jendt
1000 DO i=istrt,iendt
1001 h(i,j)=50.0_r8
1002 END DO
1003 END DO
1004#elif defined OVERFLOW
1005 val1=200.0_r8
1006 DO j=jstrt,jendt
1007 DO i=istrt,iendt
1008 h(i,j)=val1+0.5_r8*(depth-val1)* &
1009 & (1.0_r8+tanh((yr(i,j)-100000.0_r8)/20000.0_r8))
1010 END DO
1011 END DO
1012#elif defined RIVERPLUME1
1013 DO j=jstrt,jendt
1014 DO i=istrt,min(5,iendt)
1015 h(i,j)=15.0_r8
1016 END DO
1017 DO i=max(6,istrt),iendt
1018 h(i,j)=depth+real(lm(ng)-i,r8)*(15.0_r8-depth)/ &
1019 & real(lm(ng)-6,r8)
1020 END DO
1021 END DO
1022#elif defined RIVERPLUME2
1023 DO j=jstrt,jendt
1024 DO i=istrt,min(5,iendt)
1025 h(i,j)=15.0_r8
1026 END DO
1027 DO i=max(6,istrt),iendt
1028 h(i,j)=depth+real(lm(ng)-i,r8)*(15.0_r8-depth)/ &
1029 & real(lm(ng)-6,r8)
1030 END DO
1031 END DO
1032#elif defined SEAMOUNT
1033 DO j=jstrt,jendt
1034 DO i=istrt,iendt
1035 val1=(xr(i,j)-0.5_r8*xsize)/40000.0_r8
1036 val2=(yr(i,j)-0.5_r8*esize)/40000.0_r8
1037 h(i,j)=depth-4500.0_r8*exp(-(val1*val1+val2*val2))
1038 END DO
1039 END DO
1040#elif defined SED_TOY
1041 DO j=jstrt,jendt
1042 DO i=istrt,iendt
1043 h(i,j)=20.0_r8
1044 END DO
1045 END DO
1046#elif defined SHOREFACE
1047 DO j=jstrt,jendt
1048 DO i=istrt,iendt
1049 h(i,j)=11.75_r8-0.0125_r8*xsize/real(lm(ng)+1,r8)*real(i,r8)
1050 END DO
1051 END DO
1052#elif defined TEST_CHAN
1053 DO j=jstrt,jendt
1054 DO i=istrt,iendt
1055 h(i,j)=10.0_r8+0.4040_r8*real(i,r8)/real(lm(ng)+1,r8)
1056 END DO
1057 END DO
1058#elif defined UPWELLING
1059 IF (nsperiodic(ng)) THEN
1060 DO i=istrt,iendt
1061 IF (i.le.lm(ng)/2) THEN
1062 val1=real(i,r8)
1063 ELSE
1064 val1=real(lm(ng)+1-i,r8)
1065 END IF
1066 val2=min(depth,84.5_r8+66.526_r8*tanh((val1-10.0_r8)/7.0_r8))
1067 DO j=jstrt,jendt
1068 h(i,j)=val2
1069 END DO
1070 END DO
1071 ELSE IF (ewperiodic(ng)) THEN
1072 DO j=jstrt,jendt
1073 IF (j.le.mm(ng)/2) THEN
1074 val1=real(j,r8)
1075 ELSE
1076 val1=real(mm(ng)+1-j,r8)
1077 END IF
1078 val2=min(depth,84.5_r8+66.526_r8*tanh((val1-10.0_r8)/7.0_r8))
1079 DO i=istrt,iendt
1080 h(i,j)=val2
1081 END DO
1082 END DO
1083 END IF
1084#elif defined WEDDELL
1085 val1=98.80_r8
1086 val2=0.8270_r8
1087 DO k=-1,26
1088 xwrk(k)=real(k-1,r8)*15.0_r8*1000.0_r8
1089 hwrk(k)=375.0_r8
1090 END DO
1091 DO k=27,232
1092 zwrk=-2.0_r8+real(k-1,r8)*0.020_r8
1093 xwrk(k)=(520.0_r8+val1+zwrk*val1+ &
1094 & val1*val2*log(cosh(zwrk)))*1000.0_r8
1095 hwrk(k)=-75.0_r8+2198.0_r8*(1.0_r8+val2*tanh(zwrk))
1096 END DO
1097 DO k=233,235
1098 xwrk(k)=(850.0_r8+real(k-228,r8)*50.0_r8)*1000.0_r8
1099 hwrk(k)=4000.0_r8
1100 END DO
1101 DO j=jstrt,jendt
1102 DO i=istrt,iendt
1103 h(i,j)=375.0_r8
1104 DO k=1,234
1105 IF ((xwrk(k).le.xr(i,1)).and.(xr(i,1).lt.xwrk(k+1))) THEN
1106 cff=1.0_r8/(xwrk(k+1)-xwrk(k))
1107 h(i,j)=cff*(xwrk(k+1)-xr(i,j))*hwrk(k )+ &
1108 & cff*(xr(i,j)-xwrk(k ))*hwrk(k+1)
1109 END IF
1110 END DO
1111 END DO
1112 END DO
1113#elif defined WINDBASIN
1114 DO i=istrt,iendt
1115 ival=int(0.03_r8*real(lm(ng)+1,r8))
1116 IF (i.lt.ival) THEN
1117 val1=1.0_r8-(real((i+1)-ival,r8)/real(ival,r8))**2
1118 ELSE IF ((lm(ng)+1-i).lt.ival) THEN
1119 val1=1.0_r8-(real((lm(ng)+1-i)-ival,r8)/real(ival,r8))**2
1120 ELSE
1121 val1=1.0_r8
1122 END IF
1123 DO j=jstrt,jendt
1124 val2=2.0_r8*real(j-(mm(ng)+1)/2,r8)/real(mm(ng)+1,r8)
1125 h(i,j)=depth*(0.08_r8+0.92_r8*val1*(1.0_r8-val2*val2))
1126 END DO
1127 END DO
1128#else
1129 DO j=jstrt,jendt
1130 DO i=istrt,iendt
1131 h(i,j)=depth
1132 END DO
1133 END DO
1134#endif
1135!
1136! Report Statistics.
1137!
1138 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(15), 0, &
1139 & lbi, ubi, lbj, ubj, h)
1140 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1141 WRITE (stdout,10) 'bathymetry at RHO-points: h', &
1142 & ng, stats(15)%min, stats(15)%max
1143 END IF
1144 hmin(ng)=stats(15)%min
1145 hmax(ng)=stats(15)%max
1146!
1147! Exchange boundary data.
1148!
1149 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1150 CALL exchange_r2d_tile (ng, tile, &
1151 & lbi, ubi, lbj, ubj, &
1152 & h)
1153 END IF
1154
1155#ifdef DISTRIBUTE
1156 CALL mp_exchange2d (ng, tile, model, 1, &
1157 & lbi, ubi, lbj, ubj, &
1158 & nghostpoints, &
1159 & ewperiodic(ng), nsperiodic(ng), &
1160 & h)
1161#endif
1162#ifdef ICESHELF
1163!
1164!-----------------------------------------------------------------------
1165! Set depth of ice shelf (meters; negative) at RHO-points.
1166!-----------------------------------------------------------------------
1167!
1168# ifdef WEDDELL
1169 val1=340.0_r8
1170 val2=val1/16.0_r8
1171 DO j=jstrt,jendt
1172 DO i=istrt,iendt
1173 IF (i.gt.20) THEN
1174 zice(i,j)=0.0_r8
1175 ELSE IF (i.gt.4) THEN
1176 zice(i,j)=-val1+real(i-1,r8)*val2
1177 ELSE
1178 zice(i,j)=-val1
1179 END IF
1180 END DO
1181 END DO
1182# else
1183 DO j=jstrt,jendt
1184 DO i=istrt,iendt
1185 zice(i,j)=0.0_r8
1186 END DO
1187 END DO
1188# endif
1189!
1190! Report Statistics.
1191!
1192 CALL stats_2dfld (ng, tile, inlm, r2dvar, stats(16), 0, &
1193 & lbi, ubi, lbj, ubj, zice)
1194 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1195 WRITE (stdout,10) 'ice shelf thickness: zice', &
1196 & ng, stats(16)%min, stats(16)%max
1197 END IF
1198!
1199! Exchange boundary data.
1200!
1201 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1202 CALL exchange_r2d_tile (ng, tile, &
1203 & lbi, ubi, lbj, ubj, &
1204 & zice)
1205 END IF
1206
1207# ifdef DISTRIBUTE
1208 CALL mp_exchange2d (ng, tile, model, 1, &
1209 & lbi, ubi, lbj, ubj, &
1210 & nghostpoints, &
1211 & ewperiodic(ng), nsperiodic(ng), &
1212 & zice)
1213# endif
1214#endif
1215!
1216 10 FORMAT (3x,' ANA_GRID - ',a,/,19x, &
1217 & '(Grid = ',i2.2,', Min = ',1p,e15.8,0p, &
1218 & ' Max = ',1p,e15.8,0p,')')
1219!
1220 RETURN
1221 END SUBROUTINE ana_grid_tile
subroutine ana_grid_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, angler, dmde, dndx, zice, lonp, lonr, lonu, lonv, latp, latr, latu, latv, xp, xr, xu, xv, yp, yr, yu, yv, pn, pm, f, h)
Definition ana_grid.h:119
subroutine ana_grid(ng, tile, model)
Definition ana_grid.h:3
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer stdout
logical lanafile
character(len=256), dimension(39) ananame
integer, parameter inlm
Definition mod_param.F:662
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter p2dvar
Definition mod_param.F:716
integer, dimension(:), allocatable mm
Definition mod_param.F:456
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
logical spherical
real(dp), dimension(:), allocatable hmin
real(r8), dimension(:), allocatable el
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), parameter large
real(dp) eradius
real(dp), parameter deg2rad
real(dp) g
real(r8), dimension(:), allocatable xl
real(dp), dimension(:), allocatable hmax
real(dp), parameter pi
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public stats_2dfld(ng, tile, model, gtype, s, extract_flag, lbi, ubi, lbj, ubj, f, fmask, debug)
Definition stats.F:47