ROMS
Loading...
Searching...
No Matches
stiffness.F
Go to the documentation of this file.
1#include "cppdefs.h"
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 Hernan G. Arango !
9!========================================== Alexander F. Shchepetkin ===
10! !
11! This routine surveys the 3D grid in order to determine maximum !
12! grid stiffness ratio: !
13! !
14! z(i,j,k)-z(i-1,j,k)+z(i,j,k-1)-z(i-1,j,k-1) !
15! r_x = --------------------------------------------- !
16! z(i,j,k)+z(i-1,j,k)-z(i,j,k-1)-z(i-1,j,k-1) !
17! !
18! This is done for diagnostic purposes and it does not affect the !
19! computations. !
20! !
21!=======================================================================
22!
23 implicit none
24
25 PRIVATE
26 PUBLIC :: stiffness
27
28 CONTAINS
29!
30!***********************************************************************
31 SUBROUTINE stiffness (ng, tile, model)
32!***********************************************************************
33!
34 USE mod_param
35 USE mod_grid
36 USE mod_ocean
37!
38! Imported variable declarations.
39!
40 integer, intent(in) :: ng, tile, model
41!
42! Local variable declarations.
43!
44#include "tile.h"
45!
46 CALL stiffness_tile (ng, tile, model, &
47 & lbi, ubi, lbj, ubj, &
48 & imins, imaxs, jmins, jmaxs, &
49#ifdef MASKING
50 & grid(ng) % rmask, &
51 & grid(ng) % umask, &
52 & grid(ng) % vmask, &
53#endif
54 & grid(ng) % h, &
55 & grid(ng) % omn, &
56#ifdef SOLVE3D
57 & grid(ng) % Hz, &
58 & grid(ng) % z_w, &
59#endif
60 & ocean(ng)% zeta)
61 RETURN
62 END SUBROUTINE stiffness
63!
64!***********************************************************************
65 SUBROUTINE stiffness_tile (ng, tile, model, &
66 & LBi, UBi, LBj, UBj, &
67 & IminS, ImaxS, JminS, JmaxS, &
68#ifdef MASKING
69 & rmask, umask, vmask, &
70#endif
71 & h, omn, &
72#ifdef SOLVE3D
73 & Hz, z_w, &
74#endif
75 & zeta)
76!***********************************************************************
77!
78 USE mod_param
79 USE mod_parallel
80 USE mod_iounits
81 USE mod_scalars
82#ifdef DISTRIBUTE
83!
84 USE distribute_mod, ONLY : mp_reduce
85#endif
86!
87 implicit none
88!
89! Imported variable declarations.
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! Local variable declarations.
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! Compute grid stiffness.
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))))
162 DO k=1,n(ng)
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))))
180 DO k=1,n(ng)
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! Compute initial basin volume and grid cell minimum and maximum volumes.
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
202 DO k=1,n(ng)
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! Compute global values.
237!-----------------------------------------------------------------------
238!
239#ifdef DISTRIBUTE
240 nsub=1 ! distributed-memory
241#else
242 IF (domain(ng)%SouthWest_Corner(tile).and. &
243 & domain(ng)%NorthEast_Corner(tile)) THEN
244 nsub=1 ! non-tiled application
245 ELSE
246 nsub=ntilex(ng)*ntilee(ng) ! tiled application
247 END IF
248#endif
249!$OMP CRITICAL (R_FACTOR)
250 totvolume(ng)=totvolume(ng)+my_volume0
251 minvolume(ng)=min(minvolume(ng),my_volume1)
252 maxvolume(ng)=max(maxvolume(ng),my_volume2)
253#ifdef SOLVE3D
254 rx0(ng)=max(rx0(ng),my_rx0)
255 rx1(ng)=max(rx1(ng),my_rx1)
256#endif
258 IF (tile_count.eq.nsub) THEN
259 tile_count=0
260#ifdef DISTRIBUTE
261 rbuffer(1)=totvolume(ng)
262 rbuffer(2)=minvolume(ng)
263 rbuffer(3)=maxvolume(ng)
264# ifdef SOLVE3D
265 rbuffer(4)=rx0(ng)
266 rbuffer(5)=rx1(ng)
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
280 totvolume(ng)=rbuffer(1)
281 minvolume(ng)=rbuffer(2)
282 maxvolume(ng)=rbuffer(3)
283# ifdef SOLVE3D
284 rx0(ng)=rbuffer(4)
285 rx1(ng)=rbuffer(5)
286# endif
287#endif
288 IF (master.and.lwrtinfo(ng)) THEN
289 WRITE (stdout,10) ng
290 10 FORMAT (/,' Basin information for Grid ',i2.2,':',/)
291#ifdef SOLVE3D
292 WRITE (stdout,20) rx0(ng), rx1(ng)
293 20 FORMAT (' Maximum grid stiffness ratios: rx0 = ',1pe14.6, &
294 & ' (Beckmann and Haidvogel)',/,t34,'rx1 = ',1pe14.6, &
295 & ' (Haney)')
296#endif
297 IF (minvolume(ng).ne.0.0_r8) THEN
298 ratio=maxvolume(ng)/minvolume(ng)
299 ELSE
300 ratio=0.0_r8
301 END IF
302 WRITE (stdout,30) totvolume(ng), minvolume(ng), &
303 & maxvolume(ng), ratio
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!$OMP END CRITICAL (R_FACTOR)
311 RETURN
312 END SUBROUTINE stiffness_tile
313
314 END MODULE stiffness_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer stdout
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
logical master
integer tile_count
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
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
subroutine, public stiffness(ng, tile, model)
Definition stiffness.F:32
subroutine stiffness_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, rmask, umask, vmask, h, omn, hz, z_w, zeta)
Definition stiffness.F:76