81
82
87
88# ifdef DISTRIBUTE
89
91# endif
92
93 implicit none
94
95
96
97 integer, intent(in) :: ng, tile
98 integer, intent(in) :: LBi, UBi, LBj, UBj
99 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
100 integer, intent(in) :: nstp, kstp
101
102# ifdef ASSUMED_SHAPE
103 real(r8), intent(in) :: h(LBi:,LBj:)
104 real(r8), intent(in) :: omn(LBi:,LBj:)
105# ifdef SOLVE3D
106 real(r8), intent(in) :: ad_Hz(LBi:,LBj:,:)
107 real(r8), intent(in) :: ad_z_r(LBi:,LBj:,:)
108 real(r8), intent(in) :: ad_z_w(LBi:,LBj:,0:)
109 real(r8), intent(in) :: ad_rho(LBi:,LBj:,:)
110 real(r8), intent(in) :: ad_u(LBi:,LBj:,:,:)
111 real(r8), intent(in) :: ad_v(LBi:,LBj:,:,:)
112# endif
113 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
114 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
115 real(r8), intent(in) :: ad_zeta(LBi:,LBj:,:)
116# else
117 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
118 real(r8), intent(in) :: omn(LBi:UBi,LBj:UBj)
119# ifdef SOLVE3D
120 real(r8), intent(in) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
121 real(r8), intent(in) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
122 real(r8), intent(in) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
123 real(r8), intent(in) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
124 real(r8), intent(in) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
125 real(r8), intent(in) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
126# endif
127 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,:)
128 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,:)
129 real(r8), intent(in) :: ad_zeta(LBi:UBi,LBj:UBj,:)
130# endif
131
132
133
134 integer :: NSUB, i, j, k, trd
135 integer :: idia, istep
136
137# ifdef DISTRIBUTE
138 real(r8), dimension(3) :: rbuffer
139 character (len=3), dimension(3) :: op_handle
140# else
141 integer :: my_threadnum
142# endif
143
144 real(r8) :: cff, my_avgke, my_avgpe, my_volume
145
146 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ke2d
147 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: pe2d
148
149 character (len=8 ) :: kechar, pechar
150 character (len=22) :: DateTime
151
152# include "set_bounds.h"
153
154
155
156 ke2d(imins:imaxs,jmins:jmaxs)=0.0_r8
157 pe2d(imins:imaxs,jmins:jmaxs)=0.0_r8
158
159
160
161
162
163
164
165
166
167
168 istep=int(mod(real(
iic(ng)-1,r8),1.0e+10_r8))
169# ifdef SOLVE3D
170 idia=nstp
171# else
172 idia=kstp
173# endif
175
176
177
178 IF (mod(istep,
ninfo(ng)).eq.0)
THEN
179 DO j=jstr,jend
180# ifdef SOLVE3D
181 DO i=istr,iend
182 pe2d(i,j)=0.5_r8*
g*ad_z_w(i,j,
n(ng))*ad_z_w(i,j,
n(ng))
183 END DO
186 DO i=istr,iend
187 ke2d(i,j)=ke2d(i,j)+ &
188 & ad_hz(i,j,k)* &
189 & 0.25_r8*(ad_u(i ,j,k,idia)*ad_u(i ,j,k,idia)+ &
190 & ad_u(i+1,j,k,idia)*ad_u(i+1,j,k,idia)+ &
191 & ad_v(i,j ,k,idia)*ad_v(i,j ,k,idia)+ &
192 & ad_v(i,j+1,k,idia)*ad_v(i,j+1,k,idia))
193 pe2d(i,j)=pe2d(i,j)+ &
194 & cff*ad_hz(i,j,k)*(ad_rho(i,j,k)+1000.0_r8)* &
195 & (ad_z_r(i,j,k)-ad_z_w(i,j,0))
196 END DO
197 END DO
198# else
200 DO i=istr,iend
201 ke2d(i,j)=(ad_zeta(i,j,idia)+h(i,j))* &
202 & 0.25_r8*(ad_ubar(i ,j,idia)*ad_ubar(i ,j,idia)+ &
203 & ad_ubar(i+1,j,idia)*ad_ubar(i+1,j,idia)+ &
204 & ad_vbar(i,j ,idia)*ad_vbar(i,j ,idia)+ &
205 & ad_vbar(i,j+1,idia)*ad_vbar(i,j+1,idia))
206 pe2d(i,j)=cff*ad_zeta(i,j,idia)*ad_zeta(i,j,idia)
207 END DO
208# endif
209 END DO
210
211
212
213
214
215
216
217
218
219 DO j=jstr,jend
220 DO i=istr,iend
221# ifdef SOLVE3D
222
223
224# else
225
226
227# endif
228 pe2d(i,jend+1)=pe2d(i,jend+1)+omn(i,j)*h(i,j)
229 pe2d(i,jstr-1)=pe2d(i,jstr-1)+omn(i,j)*pe2d(i,j)
230 ke2d(i,jstr-1)=ke2d(i,jstr-1)+omn(i,j)*ke2d(i,j)
231 END DO
232 END DO
233 my_volume=0.0_r8
234 my_avgpe=0.0_r8
235 my_avgke=0.0_r8
236 DO i=istr,iend
237 my_volume=my_volume+pe2d(i,jend+1)
238 my_avgpe =my_avgpe +pe2d(i,jstr-1)
239 my_avgke =my_avgke +ke2d(i,jstr-1)
240 END DO
241
242
243
244
245
246
247
248# ifdef DISTRIBUTE
249 nsub=1
250# else
251 IF (
domain(ng)%SouthWest_Corner(tile).and. &
252 &
domain(ng)%NorthEast_Corner(tile))
THEN
253 nsub=1
254 ELSE
256 END IF
257# endif
258
263 END IF
270# ifdef DISTRIBUTE
274 op_handle(1)='SUM'
275 op_handle(2)='SUM'
276 op_handle(3)='SUM'
282# else
284# endif
288
289
290
293 IF (
master.and.(ng.eq.1))
THEN
294 WRITE (
stdout,10)
'TIME-STEP',
'YYYY-MM-DD hh:mm:ss.ss', &
295 & 'KINETIC_ENRG', 'POTEN_ENRG', &
296# ifdef NESTING
297 & 'TOTAL_ENRG', 'NET_VOLUME', 'Grid'
298# else
299 & 'TOTAL_ENRG', 'NET_VOLUME'
300# endif
301# ifdef NESTING
302 10 FORMAT (/,1x,a,1x,a,2x,a,3x,a,4x,a,4x,a,2x,a)
303# else
304 10 FORMAT (/,1x,a,1x,a,2x,a,3x,a,4x,a,4x,a)
305# endif
306 END IF
307 END IF
308
310 WRITE (
stdout,20) istep, datetime, &
311# ifdef NESTING
313# else
315# endif
316# ifdef NESTING
317 20 FORMAT (i10,1x,a,4(1pe14.6),2x,i2.2)
318# else
319 20 FORMAT (i10,1x,a,4(1pe14.6))
320# endif
322 END IF
323
324
325
326 WRITE (kechar,
'(1pe8.1)')
avgke
327 WRITE (pechar,
'(1pe8.1)')
avgpe
328 DO i=1,8
329 IF ((kechar(i:i).eq.'N').or.(pechar(i:i).eq.'N').or. &
330 & (kechar(i:i).eq.'n').or.(pechar(i:i).eq.'n').or. &
331 & (kechar(i:i).eq.'*').or.(pechar(i:i).eq.'*')) THEN
333 END IF
334 END DO
335 END IF
336
337 END IF
338
339 RETURN
integer function my_threadnum()
integer, dimension(:), allocatable n
integer, dimension(:), allocatable ntilex
type(t_domain), dimension(:), allocatable domain
integer, dimension(:), allocatable ntilee
integer, dimension(:), allocatable ninfo
integer, dimension(:), allocatable iic
integer, dimension(:), allocatable first_time
character(len=22), dimension(:), allocatable time_code