76
77
82#ifdef DISTRIBUTE
83
85#endif
86
87 implicit none
88
89
90
91 integer, intent(in) :: ng, tile, model
92 integer, intent(in) :: LBi, UBi, LBj, UBj
93 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
94
95#ifdef ASSUMED_SHAPE
96# ifdef MASKING
97 real(r8), intent(in) :: rmask(LBi:,LBj:)
98 real(r8), intent(in) :: umask(LBi:,LBj:)
99 real(r8), intent(in) :: vmask(LBi:,LBj:)
100# endif
101 real(r8), intent(in) :: h(LBi:,LBj:)
102 real(r8), intent(in) :: omn(LBi:,LBj:)
103# ifdef SOLVE3D
104 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
105 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
106# endif
107 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
108#else
109# ifdef MASKING
110 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
111 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
112 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
113# endif
114 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
115 real(r8), intent(in) :: omn(LBi:UBi,LBj:UBj)
116# ifdef SOLVE3D
117 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
118 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
119# endif
120 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
121#endif
122
123
124
125 integer :: NSUB, i, j, k
126
127 real(r8) :: cff, ratio
128
129#ifdef SOLVE3D
130 real(r8) :: my_rx0, my_rx1
131#endif
132 real(r8) :: my_volume0, my_volume1, my_volume2
133
134#ifdef DISTRIBUTE
135# ifdef SOLVE3D
136 real(r8), dimension(5) :: rbuffer
137 character (len=3), dimension(5) :: op_handle
138# else
139 real(r8), dimension(3) :: rbuffer
140 character (len=3), dimension(3) :: op_handle
141# endif
142#endif
143
144#include "set_bounds.h"
145
146#ifdef SOLVE3D
147
148
149
150
151
152 my_rx0=0.0_r8
153 my_rx1=0.0_r8
154
155 DO j=jstr,jend
156 DO i=istru,iend
157# ifdef MASKING
158 IF (umask(i,j).gt.0.0_r8) THEN
159# endif
160 my_rx0=max(my_rx0,abs((z_w(i,j,0)-z_w(i-1,j,0))/ &
161 & (z_w(i,j,0)+z_w(i-1,j,0))))
163 my_rx1=max(my_rx1,abs((z_w(i,j,k )-z_w(i-1,j,k )+ &
164 & z_w(i,j,k-1)-z_w(i-1,j,k-1))/ &
165 & (z_w(i,j,k )+z_w(i-1,j,k )- &
166 & z_w(i,j,k-1)-z_w(i-1,j,k-1))))
167 END DO
168# ifdef MASKING
169 END IF
170# endif
171 END DO
172 END DO
173 DO j=jstrv,jend
174 DO i=istr,iend
175# ifdef MASKING
176 IF (vmask(i,j).gt.0.0_r8) THEN
177# endif
178 my_rx0=max(my_rx0,abs((z_w(i,j,0)-z_w(i,j-1,0))/ &
179 & (z_w(i,j,0)+z_w(i,j-1,0))))
181 my_rx1=max(my_rx1,abs((z_w(i,j,k )-z_w(i,j-1,k )+ &
182 & z_w(i,j,k-1)-z_w(i,j-1,k-1))/ &
183 & (z_w(i,j,k )+z_w(i,j-1,k )- &
184 & z_w(i,j,k-1)-z_w(i,j-1,k-1))))
185 END DO
186# ifdef MASKING
187 END IF
188# endif
189 END DO
190 END DO
191#endif
192
193
194
195
196
197 my_volume0=0.0_r8
198 my_volume1=1.0e+20_r8
199 my_volume2=0.0_r8
200
201#ifdef SOLVE3D
203 DO j=jstr,jend
204 DO i=istr,iend
205# ifdef MASKING
206 IF (rmask(i,j).gt.0.0_r8) THEN
207# endif
208 cff=omn(i,j)*hz(i,j,k)
209 my_volume0=my_volume0+cff
210 my_volume1=min(my_volume1,cff)
211 my_volume2=max(my_volume2,cff)
212# ifdef MASKING
213 END IF
214# endif
215 END DO
216 END DO
217 END DO
218#else
219 DO j=jstr,jend
220 DO i=istr,iend
221# ifdef MASKING
222 IF (rmask(i,j).gt.0.0_r8) THEN
223# endif
224 cff=omn(i,j)*(zeta(i,j,1)+h(i,j))
225 my_volume0=my_volume0+cff
226 my_volume1=min(my_volume1,cff)
227 my_volume2=max(my_volume2,cff)
228# ifdef MASKING
229 END IF
230# endif
231 END DO
232 END DO
233#endif
234
235
236
237
238
239#ifdef DISTRIBUTE
240 nsub=1
241#else
242 IF (
domain(ng)%SouthWest_Corner(tile).and. &
243 &
domain(ng)%NorthEast_Corner(tile))
THEN
244 nsub=1
245 ELSE
247 END IF
248#endif
249
253#ifdef SOLVE3D
254 rx0(ng)=max(
rx0(ng),my_rx0)
255 rx1(ng)=max(
rx1(ng),my_rx1)
256#endif
260#ifdef DISTRIBUTE
264# ifdef SOLVE3D
267# endif
268 op_handle(1)='SUM'
269 op_handle(2)='MIN'
270 op_handle(3)='MAX'
271# ifdef SOLVE3D
272 op_handle(4)='MAX'
273 op_handle(5)='MAX'
274# endif
275# ifdef SOLVE3D
276 CALL mp_reduce (ng, model, 5, rbuffer, op_handle)
277# else
278 CALL mp_reduce (ng, model, 3, rbuffer, op_handle)
279# endif
283# ifdef SOLVE3D
286# endif
287#endif
290 10 FORMAT (/,' Basin information for Grid ',i2.2,':',/)
291#ifdef SOLVE3D
293 20 FORMAT (' Maximum grid stiffness ratios: rx0 = ',1pe14.6, &
294 & ' (Beckmann and Haidvogel)',/,t34,'rx1 = ',1pe14.6, &
295 & ' (Haney)')
296#endif
299 ELSE
300 ratio=0.0_r8
301 END IF
304 30 FORMAT (/,' Initial domain volumes: TotVolume = ',1p,e17.10, &
305 & 0p,' m3',/,t26,'MinCellVol = ',1p,e17.10,0p,' m3', &
306 & /,t26, 'MaxCellVol = ',1p,e17.10,0p,' m3', &
307 & /,t26, ' Max/Min = ',1p,e17.10,0p,/)
308 END IF
309 END IF
310
311 RETURN
integer, dimension(:), allocatable n
integer, dimension(:), allocatable ntilex
type(t_domain), dimension(:), allocatable domain
integer, dimension(:), allocatable ntilee
real(dp), dimension(:), allocatable totvolume
real(dp), dimension(:), allocatable minvolume
real(dp), dimension(:), allocatable maxvolume
logical, dimension(:), allocatable lwrtinfo
real(dp), dimension(:), allocatable rx1
real(dp), dimension(:), allocatable rx0