ROMS
Loading...
Searching...
No Matches
wpoints.F File Reference
#include "cppdefs.h"
#include "tile.h"
#include "set_bounds.h"
Include dependency graph for wpoints.F:

Go to the source code of this file.

Macros

#define IR_RANGE   IstrT,IendT
 
#define IU_RANGE   IstrP,IendT
 
#define JR_RANGE   JstrT,JendT
 
#define JV_RANGE   JstrP,JendT
 

Functions/Subroutines

program __wpoints_f__
 
subroutine wpoints (ng, tile, model)
 
subroutine wpoints_tile (ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, pmask_full, rmask_full, umask_full, vmask_full, ijwaterr, ijwateru, ijwaterv, h)
 

Macro Definition Documentation

◆ IR_RANGE

#define IR_RANGE   IstrT,IendT

◆ IU_RANGE

#define IU_RANGE   IstrP,IendT

◆ JR_RANGE

#define JR_RANGE   JstrT,JendT

◆ JV_RANGE

#define JV_RANGE   JstrP,JendT

Function/Subroutine Documentation

◆ __wpoints_f__()

program __wpoints_f__

Definition at line 3 of file wpoints.F.

References wpoints().

Here is the call graph for this function:

◆ wpoints()

subroutine __wpoints_f__::wpoints ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model )
private

Definition at line 39 of file wpoints.F.

40!***********************************************************************
41!
42 USE mod_param
43 USE mod_grid
44!
45 implicit none
46!
47! Imported variable declarations.
48!
49 integer, intent(in) :: ng, tile, model
50!
51! Local variable declarations.
52!
53# include "tile.h"
54!
55 CALL wpoints_tile (ng, tile, model, &
56 & lbi, ubi, lbj, ubj, &
57 & imins, imaxs, jmins, jmaxs, &
58# ifdef MASKING
59 & grid(ng) % pmask_full, &
60 & grid(ng) % rmask_full, &
61 & grid(ng) % umask_full, &
62 & grid(ng) % vmask_full, &
63# ifdef PROPAGATOR
64 & grid(ng) % IJwaterR, &
65 & grid(ng) % IJwaterU, &
66 & grid(ng) % IJwaterV, &
67# endif
68# endif
69 & grid(ng) % h)
70
71 RETURN
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
subroutine wpoints_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, pmask_full, rmask_full, umask_full, vmask_full, ijwaterr, ijwateru, ijwaterv, h)
Definition wpoints.F:86

References mod_grid::grid, and wpoints_tile().

Referenced by __wpoints_f__(), ad_initial(), initial(), rp_initial(), and tl_initial().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ wpoints_tile()

subroutine __wpoints_f__::wpoints_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) model,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
real(r8), dimension(lbi:,lbj:), intent(in) pmask_full,
real(r8), dimension(lbi:,lbj:), intent(in) rmask_full,
real(r8), dimension(lbi:,lbj:), intent(in) umask_full,
real(r8), dimension(lbi:,lbj:), intent(in) vmask_full,
integer, dimension(lbi:,lbj:), intent(out) ijwaterr,
integer, dimension(lbi:,lbj:), intent(out) ijwateru,
integer, dimension(lbi:,lbj:), intent(out) ijwaterv,
real(r8), dimension(lbi:,lbj:), intent(in) h )
private

Definition at line 75 of file wpoints.F.

86!***********************************************************************
87!
88 USE mod_param
89 USE mod_parallel
90 USE mod_iounits
91 USE mod_ncparam
92# if defined PROPAGATOR && defined CHECKPOINTING && \
93 defined pio_lib && defined distribute
95# endif
96 USE mod_scalars
97# ifdef PROPAGATOR
98 USE mod_storage
99# endif
100# ifdef DISTRIBUTE
101!
103# if defined PROPAGATOR && defined CHECKPOINTING && defined PIO_LIB
104 USE set_pio_mod, ONLY : state_iodecomp
105# endif
106# endif
107!
108 implicit none
109!
110! Imported variable declarations.
111!
112 integer, intent(in) :: ng, tile, model
113 integer, intent(in) :: LBi, UBi, LBj, UBj
114 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
115!
116# ifdef ASSUMED_SHAPE
117# ifdef MASKING
118# ifdef PROPAGATOR
119 integer, intent(out) :: IJwaterR(LBi:,LBj:)
120 integer, intent(out) :: IJwaterU(LBi:,LBj:)
121 integer, intent(out) :: IJwaterV(LBi:,LBj:)
122# endif
123 real(r8), intent(in) :: rmask_full(LBi:,LBj:)
124 real(r8), intent(in) :: pmask_full(LBi:,LBj:)
125 real(r8), intent(in) :: umask_full(LBi:,LBj:)
126 real(r8), intent(in) :: vmask_full(LBi:,LBj:)
127# endif
128 real(r8), intent(in) :: h(LBi:,LBj:)
129# else
130# ifdef MASKING
131# ifdef PROPAGATOR
132 integer, intent(out) :: IJwaterR(LBi:UBi,LBj:UBj)
133 integer, intent(out) :: IJwaterU(LBi:UBi,LBj:UBj)
134 integer, intent(out) :: IJwaterV(LBi:UBi,LBj:UBj)
135# endif
136 real(r8), intent(in) :: rmask_full(LBi:UBi,LBj:UBj)
137 real(r8), intent(in) :: pmask_full(LBi:UBi,LBj:UBj)
138 real(r8), intent(in) :: umask_full(LBi:UBi,LBj:UBj)
139 real(r8), intent(in) :: vmask_full(LBi:UBi,LBj:UBj)
140# endif
141 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
142# endif
143!
144! Local variable declarations.
145!
146# if defined PROPAGATOR && defined CHECKPOINTING && \
147 defined pio_lib && defined distribute
148 logical, save :: first = .true.
149!
150# endif
151 integer :: my_Nxyp, my_Nxyr, my_Nxyu, my_Nxyv
152# ifdef PROPAGATOR
153 integer :: my_NwaterR, my_NwaterU, my_NwaterV
154# endif
155
156# ifdef PROPAGATOR
157 integer :: Imin, Imax, Jmin, Jmax
158# endif
159 integer :: Uoff, Voff
160 integer :: NSUB, Npts, i, ic, ij, j
161# ifdef DISTRIBUTE
162# ifdef PROPAGATOR
163 integer :: block_size
164
165 real(r8), dimension(3) :: wp_buffer
166 character (len=3), dimension(3) :: wp_handle
167# endif
168# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
169 real(r8), dimension(4) :: io_buffer
170 character (len=3), dimension(4) :: io_handle
171# endif
172# endif
173# if defined READ_WATER || defined PROPAGATOR
174 real(r8), dimension((Lm(ng)+2)*(Mm(ng)+2)) :: mask
175# endif
176
177# include "set_bounds.h"
178
179# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
180!
181!-----------------------------------------------------------------------
182! Determine number of water points.
183!-----------------------------------------------------------------------
184!
185! Determine interior points for IO purposes.
186!
187 my_nxyr=0
188 my_nxyp=0
189 my_nxyu=0
190 my_nxyv=0
191
192 DO j=jstr,jend
193 DO i=istr,iend
194 IF (pmask_full(i,j).gt.0.0_r8) my_nxyp=my_nxyp+1
195 IF (rmask_full(i,j).gt.0.0_r8) my_nxyr=my_nxyr+1
196 IF (umask_full(i,j).gt.0.0_r8) my_nxyu=my_nxyu+1
197 IF (vmask_full(i,j).gt.0.0_r8) my_nxyv=my_nxyv+1
198 END DO
199 END DO
200!
201! Determine boundary points for IO purposes. The assigment
202! below is done to account for periodicity.
203!
204 IF (domain(ng)%Western_Edge(tile)) THEN
205 DO j=jstr,jend
206 IF (rmask_full(istr-1,j).gt.0.0_r8) my_nxyr=my_nxyr+1
207 IF (vmask_full(istr-1,j).gt.0.0_r8) my_nxyv=my_nxyv+1
208 END DO
209 END IF
210 IF (domain(ng)%Eastern_Edge(tile)) THEN
211 DO j=jstr,jend
212 IF (pmask_full(iend+1,j).gt.0.0_r8) my_nxyp=my_nxyp+1
213 IF (rmask_full(iend+1,j).gt.0.0_r8) my_nxyr=my_nxyr+1
214 IF (umask_full(iend+1,j).gt.0.0_r8) my_nxyu=my_nxyu+1
215 IF (vmask_full(iend+1,j).gt.0.0_r8) my_nxyv=my_nxyv+1
216 END DO
217 END IF
218 IF (domain(ng)%Southern_Edge(tile)) THEN
219 DO i=istr,iend
220 IF (rmask_full(i,jstr-1).gt.0.0_r8) my_nxyr=my_nxyr+1
221 IF (umask_full(i,jstr-1).gt.0.0_r8) my_nxyu=my_nxyu+1
222 END DO
223 END IF
224 IF (domain(ng)%Northern_Edge(tile)) THEN
225 DO i=istr,iend
226 IF (pmask_full(i,jend+1).gt.0.0_r8) my_nxyp=my_nxyp+1
227 IF (rmask_full(i,jend+1).gt.0.0_r8) my_nxyr=my_nxyr+1
228 IF (umask_full(i,jend+1).gt.0.0_r8) my_nxyu=my_nxyu+1
229 IF (vmask_full(i,jend+1).gt.0.0_r8) my_nxyv=my_nxyv+1
230 END DO
231 END IF
232 IF (domain(ng)%SouthWest_Corner(tile)) THEN
233 IF (rmask_full(istr-1,jstr-1).gt.0.0_r8) my_nxyr=my_nxyr+1
234 END IF
235 IF (domain(ng)%SouthEast_Corner(tile)) THEN
236 IF (rmask_full(iend+1,jstr-1).gt.0.0_r8) my_nxyr=my_nxyr+1
237 IF (umask_full(iend+1,jstr-1).gt.0.0_r8) my_nxyu=my_nxyu+1
238 END IF
239 IF (domain(ng)%NorthWest_Corner(tile)) THEN
240 IF (rmask_full(istr-1,jend+1).gt.0.0_r8) my_nxyr=my_nxyr+1
241 IF (vmask_full(istr-1,jend+1).gt.0.0_r8) my_nxyv=my_nxyv+1
242 END IF
243 IF (domain(ng)%NorthEast_Corner(tile)) THEN
244 IF (pmask_full(iend+1,jend+1).gt.0.0_r8) my_nxyp=my_nxyp+1
245 IF (rmask_full(iend+1,jend+1).gt.0.0_r8) my_nxyr=my_nxyr+1
246 IF (umask_full(iend+1,jend+1).gt.0.0_r8) my_nxyu=my_nxyu+1
247 IF (vmask_full(iend+1,jend+1).gt.0.0_r8) my_nxyv=my_nxyv+1
248 END IF
249# endif
250# ifdef PROPAGATOR
251!
252!-----------------------------------------------------------------------
253! Determine number of water points for propagator state vector.
254!-----------------------------------------------------------------------
255!
256 my_nwaterr=0
257 my_nwateru=0
258 my_nwaterv=0
259
260 DO j=jr_range
261 DO i=ir_range
262# ifdef MASKING
263 IF (rmask_full(i,j).gt.0.0_r8) THEN
264 my_nwaterr=my_nwaterr+1
265 END IF
266# else
267 my_nwaterr=my_nwaterr+1
268# endif
269 END DO
270 DO i=iu_range
271# ifdef MASKING
272 IF (umask_full(i,j).gt.0.0_r8) THEN
273 my_nwateru=my_nwateru+1
274 END IF
275# else
276 my_nwateru=my_nwateru+1
277# endif
278 END DO
279 END DO
280 DO j=jv_range
281 DO i=ir_range
282# ifdef MASKING
283 IF (vmask_full(i,j).gt.0.0_r8) THEN
284 my_nwaterv=my_nwaterv+1
285 END IF
286# else
287 my_nwaterv=my_nwaterv+1
288# endif
289 END DO
290 END DO
291# endif
292!
293! Determine global number of points (parallel reduction).
294!
295# ifdef DISTRIBUTE
296 nsub=1 ! distributed-memory
297# else
298 IF (domain(ng)%SouthWest_Corner(tile).and. &
299 & domain(ng)%NorthEast_Corner(tile)) THEN
300 nsub=1 ! non-tiled application
301 ELSE
302 nsub=ntilex(ng)*ntilee(ng) ! tiled application
303 END IF
304# endif
305!$OMP CRITICAL (MAX_WATER)
306 IF (block_count.eq.0) THEN
307# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
308 nxyp(ng)=0
309 nxyr(ng)=0
310 nxyu(ng)=0
311 nxyv(ng)=0
312# endif
313# ifdef PROPAGATOR
314 nwaterr(ng)=0
315 nwateru(ng)=0
316 nwaterv(ng)=0
317# endif
318 END IF
319# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
320 nxyp(ng)=nxyp(ng)+my_nxyp
321 nxyr(ng)=nxyr(ng)+my_nxyr
322 nxyu(ng)=nxyu(ng)+my_nxyu
323 nxyv(ng)=nxyv(ng)+my_nxyv
324# endif
325# ifdef PROPAGATOR
326 nwaterr(ng)=nwaterr(ng)+my_nwaterr
327 nwateru(ng)=nwateru(ng)+my_nwateru
328 nwaterv(ng)=nwaterv(ng)+my_nwaterv
329# endif
331 IF (block_count.eq.nsub) THEN
333# ifdef DISTRIBUTE
334# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
335 io_buffer(1)=real(nxyp(ng),r8)
336 io_buffer(2)=real(nxyr(ng),r8)
337 io_buffer(3)=real(nxyu(ng),r8)
338 io_buffer(4)=real(nxyv(ng),r8)
339 io_handle(1)='SUM'
340 io_handle(2)='SUM'
341 io_handle(3)='SUM'
342 io_handle(4)='SUM'
343 CALL mp_reduce (ng, model, 4, io_buffer, io_handle)
344 nxyp(ng)=int(io_buffer(1))
345 nxyr(ng)=int(io_buffer(2))
346 nxyu(ng)=int(io_buffer(3))
347 nxyv(ng)=int(io_buffer(4))
348# endif
349# ifdef PROPAGATOR
350 wp_buffer(1)=real(nwaterr(ng),r8)
351 wp_buffer(2)=real(nwateru(ng),r8)
352 wp_buffer(3)=real(nwaterv(ng),r8)
353 wp_handle(1)='SUM'
354 wp_handle(2)='SUM'
355 wp_handle(3)='SUM'
356 CALL mp_reduce (ng, model, 3, wp_buffer, wp_handle)
357 nwaterr(ng)=int(wp_buffer(1))
358 nwateru(ng)=int(wp_buffer(2))
359 nwaterv(ng)=int(wp_buffer(3))
360# endif
361# endif
362# if defined MASKING && (defined READ_WATER || defined WRITE_WATER)
363 iobounds(ng) % xy_psi = nxyp(ng)
364 iobounds(ng) % xy_rho = nxyr(ng)
365 iobounds(ng) % xy_u = nxyu(ng)
366 iobounds(ng) % xy_v = nxyv(ng)
367# endif
368 END IF
369!$OMP END CRITICAL (MAX_WATER)
370
371# ifdef PROPAGATOR
372!
373!-----------------------------------------------------------------------
374! Determine size of state vector.
375!-----------------------------------------------------------------------
376# if defined FORCING_SV || \
377 defined so_semi || defined stochastic_opt
378!
379! The state vector may contain free-surface, momentum, tracers,
380! surface momentum stresses, amd surface tracer fluxes.
381!
382 mstate(ng)=0
383 IF (scalars(ng)%Fstate(isfsur)) THEN
384 mstate(ng)=mstate(ng)+nwaterr(ng)
385 END IF
386# ifdef SOLVE3D
387 IF (scalars(ng)%Fstate(isuvel)) THEN
388 mstate(ng)=mstate(ng)+n(ng)*nwateru(ng)
389 END IF
390 IF (scalars(ng)%Fstate(isvvel)) THEN
391 mstate(ng)=mstate(ng)+n(ng)*nwaterv(ng)
392 END IF
393 DO i=1,nt(ng)
394 IF (scalars(ng)%Fstate(istvar(i))) THEN
395 mstate(ng)=mstate(ng)+n(ng)*nwaterr(ng)
396 END IF
397 END DO
398# else
399 IF (scalars(ng)%Fstate(isubar)) THEN
400 mstate(ng)=mstate(ng)+nwateru(ng)
401 END IF
402 IF (scalars(ng)%Fstate(isvbar)) THEN
403 mstate(ng)=mstate(ng)+nwaterv(ng)
404 END IF
405# endif
406 IF (scalars(ng)%Fstate(isustr)) THEN
407 mstate(ng)=mstate(ng)+nwateru(ng)
408 END IF
409 IF (scalars(ng)%Fstate(isvstr)) THEN
410 mstate(ng)=mstate(ng)+nwaterv(ng)
411 END IF
412# ifdef SOLVE3D
413 DO i=1,nt(ng)
414 IF (scalars(ng)%Fstate(istsur(i))) THEN
415 mstate(ng)=mstate(ng)+nwaterr(ng)
416 END IF
417 END DO
418# endif
419# else
420# ifdef SOLVE3D
421!
422! The state vector contains free-surface, 3D momentum components, and
423! tracers.
424!
425 mstate(ng)=nwaterr(ng)+ &
426 & nwateru(ng)*n(ng)+ &
427 & nwaterv(ng)*n(ng)+ &
428 & nwaterr(ng)*n(ng)*nt(ng)
429# else
430!
431! The state vector contains free-surface and 2D momentum components.
432!
433 mstate(ng)=nwaterr(ng)+ &
434 & nwateru(ng)+ &
435 & nwaterv(ng)
436# endif
437# endif
438# ifdef DISTRIBUTE
439!
440! Deternine state vector starting and ending blocking indices.
441!
443 block_size=nstate(ng)
444 IF (master) THEN
445 WRITE (stdout,10) ng, mstate(ng), nstate(ng)
446 10 FORMAT (/,' Size of FULL state arrays for propagator of grid ', &
447 & i2.2,', Mstate = ', i10,/,38x, &
448 & 'NODE partition, Nstate = ', i10,/,/, &
449 & 9x,'Node',7x,'Nstr',7x,'Nend',7x,'Size',/)
450 DO ic=0,numthreads-1
451 i=1+ic*block_size
452 j=min(mstate(ng),i+block_size-1)
453 WRITE (stdout,20) ic, i, j, j-i+1
454 20 FORMAT (11x, i2, 3(1x,i10))
455 END DO
456 WRITE (stdout,'(/)')
457 END IF
458 nstr(ng)=1+myrank*block_size
459 nend(ng)=nstr(ng)+block_size-1
460 mstate(ng)=nstate(ng)*numthreads
461# else
462 nstr(ng)=1
463 nend(ng)=mstate(ng)
464 nstate(ng)=mstate(ng)
465 IF (master) THEN
466 WRITE (stdout,10) ng, mstate(ng)
467 10 FORMAT (' Size of FULL state arrays for propagator of grid ', &
468 & i2.2, ', Mstate = ', i10,/)
469 END IF
470# endif
471# endif
472
473# if defined PIO_LIB && defined DISTRIBUTE && defined CHECKPOINTING
474!
475!-----------------------------------------------------------------------
476! Set I/O decomposition descriptors for single and double precision
477! state propagator data
478!-----------------------------------------------------------------------
479!
480 IF (first) THEN
481 IF (ng.eq.ngrids) first=.false.
482!
483 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_real, &
484 & iodesc_sp_bvec(ng), &
485 & 'Bvec', 2)
486 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_real, &
487 & iodesc_sp_resid(ng), &
488 & 'resid', 1)
489 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_real, &
490 & iodesc_sp_sworkd(ng), &
491 & 'SworkD', 1)
492!
493 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_double, &
494 & iodesc_dp_bvec(ng), &
495 & 'Bvec', 2)
496 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_double, &
497 & iodesc_dp_resid(ng), &
498 & 'resid', 1)
499 CALL state_iodecomp (ng, piosystem(ipioroms,ng), pio_double, &
500 & iodesc_dp_sworkd(ng), &
501 & 'SworkD', 1)
502 END IF
503# endif
504
505# if ((defined READ_WATER && defined DISTRIBUTE) || \
506 defined propagator) && defined masking
507!
508!-----------------------------------------------------------------------
509! If distributed-memory, set process water indices arrays.
510!-----------------------------------------------------------------------
511!
512! Set offsets for U- and V-variables due to periodic boundary
513! conditions. Recall that in East-West periodic boundary conditions
514! IstrU=1. Otherwise, IstrU=2. Similarly, in North-South periodic
515! applications IstrV=1 or else IstrV=2.
516!
517 IF (ewperiodic(ng)) THEN
518 uoff=0
519 ELSE
520 uoff=1
521 END IF
522!
523 IF (nsperiodic(ng)) THEN
524 voff=0
525 ELSE
526 voff=1
527 END IF
528
529# ifdef DISTRIBUTE
530!
531! Gather masking at RHO-points from all nodes.
532!
533 CALL mp_gather2d (ng, model, lbi, ubi, lbj, ubj, &
534 & 0, r2dvar, 1.0_r8, &
535 & rmask_full, rmask_full, npts, mask, .false.)
536 CALL mp_bcastf (ng, model, mask)
537# else
538!
539! Load masking at RHO-points in mask vector.
540!
541 ic=0
542 DO j=0,mm(ng)+1
543 DO i=0,lm(ng)+1
544 ic=ic+1
545 mask(ic)=rmask_full(i,j)
546 END DO
547 END DO
548 npts=ic
549# endif
550# if defined READ_WATER && defined DISTRIBUTE
551!
552! Determine water RHO-points for IO purposes.
553!
554 ic=0
555 DO ij=1,npts
556 IF (mask(ij).gt.0.0_r8) THEN
557 ic=ic+1
558 scalars(ng)%IJwater(ic,r2dvar)=ij
559 END IF
560 END DO
561 my_nxyr=ic
562# endif
563# ifdef PROPAGATOR
564!
565! Determine water points for the propagator state vector at RHO-points.
566!
567 ic=0
568 ij=0
569# ifdef FULL_GRID
570 imin=0
571 imax=lm(ng)+1
572 jmin=0
573 jmax=mm(ng)+1
574# else
575 imin=1
576 imax=lm(ng)
577 jmin=1
578 jmax=mm(ng)
579# endif
580 DO j=0,mm(ng)+1
581 DO i=0,lm(ng)+1
582 ij=ij+1
583 IF ((mask(ij).gt.0.0_r8).and. &
584 & (imin.le.i).and.(i.le.imax).and. &
585 & (jmin.le.j).and.(j.le.jmax)) THEN
586 ic=ic+1
587 mask(ij)=real(ic,r8)
588 ELSE
589 mask(ij)=0.0_r8
590 END IF
591 END DO
592 END DO
593!
594 ij=0
595 DO j=0,mm(ng)+1
596 DO i=0,lm(ng)+1
597 ij=ij+1
598 IF ((rilb(ng).le.i).and.(i.le.riub(ng)).and. &
599 & (rjlb(ng).le.j).and.(j.le.rjub(ng))) THEN
600 IF (mask(ij).gt.0.0_r8) THEN
601 ijwaterr(i,j)=int(mask(ij))
602 ELSE
603 ijwaterr(i,j)=0
604 END IF
605 END IF
606 END DO
607 END DO
608# endif
609
610# if defined READ_WATER && defined DISTRIBUTE
611!
612! Gather masking at PSI-points from all nodes.
613!
614 CALL mp_gather2d (ng, model, lbi, ubi, lbj, ubj, &
615 & 0, p2dvar, 1.0_r8, &
616 & pmask_full, pmask_full, npts, mask, .false.)
617 CALL mp_bcastf (ng, model, mask)
618!
619! Determine water PSI-points for IO purposes.
620!
621 ic=0
622 DO ij=1,npts
623 IF (mask(ij).gt.0.0_r8) THEN
624 ic=ic+1
625 scalars(ng)%IJwater(ic,p2dvar)=ij
626 END IF
627 END DO
628 my_nxyp=ic
629# endif
630
631# ifdef DISTRIBUTE
632!
633! Gather masking at U-points from all nodes.
634!
635 CALL mp_gather2d (ng, model, lbi, ubi, lbj, ubj, &
636 & 0, u2dvar, 1.0_r8, &
637 & umask_full, umask_full, npts, mask, .false.)
638 CALL mp_bcastf (ng, model, mask)
639# else
640!
641! Load masking at U-points in mask vector.
642!
643 ic=0
644 DO j=0,mm(ng)+1
645 DO i=1,lm(ng)+1
646 ic=ic+1
647 mask(ic)=umask_full(i,j)
648 END DO
649 END DO
650 npts=ic
651# endif
652# if defined READ_WATER && defined DISTRIBUTE
653!
654! Determine water U-points for IO purposes.
655!
656 ic=0
657 DO ij=1,npts
658 IF (mask(ij).gt.0.0_r8) THEN
659 ic=ic+1
660 scalars(ng)%IJwater(ic,u2dvar)=ij
661 END IF
662 END DO
663 my_nxyu=ic
664# endif
665# ifdef PROPAGATOR
666!
667! Determine water points for the propagator state vector at U-points.
668!
669 ic=0
670 ij=0
671# ifdef FULL_GRID
672 imin=1
673 imax=lm(ng)+1
674 jmin=0
675 jmax=mm(ng)+1
676# else
677 imin=1+uoff
678 imax=lm(ng)
679 jmin=1
680 jmax=mm(ng)
681# endif
682 DO j=0,mm(ng)+1
683 DO i=1,lm(ng)+1
684 ij=ij+1
685 IF ((mask(ij).gt.0.0_r8).and. &
686 & (imin.le.i).and.(i.le.imax).and. &
687 & (jmin.le.j).and.(j.le.jmax)) THEN
688 ic=ic+1
689 mask(ij)=real(ic,r8)
690 ELSE
691 mask(ij)=0.0_r8
692 END IF
693 END DO
694 END DO
695!
696 ij=0
697 DO j=0,mm(ng)+1
698 DO i=1,lm(ng)+1
699 ij=ij+1
700 IF ((uilb(ng).le.i).and.(i.le.uiub(ng)).and. &
701 & (ujlb(ng).le.j).and.(j.le.ujub(ng))) THEN
702 IF (mask(ij).gt.0.0_r8) THEN
703 ijwateru(i,j)=int(mask(ij))
704 ELSE
705 ijwateru(i,j)=0
706 END IF
707 END IF
708 END DO
709 END DO
710# endif
711
712# ifdef DISTRIBUTE
713!
714! Gather masking at V-points from all nodes.
715!
716 CALL mp_gather2d (ng, model, lbi, ubi, lbj, ubj, &
717 & 0, v2dvar, 1.0_r8, &
718 & vmask_full, vmask_full, npts, mask, .false.)
719 CALL mp_bcastf (ng, model, mask)
720# else
721!
722! Load masking at V-points in mask vector.
723!
724 ic=0
725 DO j=1,mm(ng)+1
726 DO i=0,lm(ng)+1
727 ic=ic+1
728 mask(ic)=vmask_full(i,j)
729 END DO
730 END DO
731 npts=ic
732# endif
733# if defined READ_WATER && defined DISTRIBUTE
734!
735! Determine water V-points for IO purposes.
736!
737 ic=0
738 DO ij=1,npts
739 IF (mask(ij).gt.0.0_r8) THEN
740 ic=ic+1
741 scalars(ng)%IJwater(ic,v2dvar)=ij
742 END IF
743 END DO
744 my_nxyv=ic
745# endif
746# ifdef PROPAGATOR
747!
748! Determine water points for the propagator state vector at V-points.
749!
750 ic=0
751 ij=0
752# ifdef FULL_GRID
753 imin=0
754 imax=lm(ng)+1
755 jmin=1
756 jmax=mm(ng)+1
757# else
758 imin=1
759 imax=lm(ng)
760 jmin=1+voff
761 jmax=mm(ng)
762# endif
763 DO j=1,mm(ng)+1
764 DO i=0,lm(ng)+1
765 ij=ij+1
766 IF ((mask(ij).gt.0.0_r8).and. &
767 & (imin.le.i).and.(i.le.imax).and. &
768 & (jmin.le.j).and.(j.le.jmax)) THEN
769 ic=ic+1
770 mask(ij)=real(ic,r8)
771 ELSE
772 mask(ij)=0.0_r8
773 END IF
774 END DO
775 END DO
776!
777 ij=0
778 DO j=1,mm(ng)+1
779 DO i=0,lm(ng)+1
780 ij=ij+1
781 IF ((vilb(ng).le.i).and.(i.le.viub(ng)).and. &
782 & (vjlb(ng).le.j).and.(j.le.vjub(ng))) THEN
783 IF (mask(ij).gt.0.0_r8) THEN
784 ijwaterv(i,j)=int(mask(ij))
785 ELSE
786 ijwaterv(i,j)=0
787 END IF
788 END IF
789 END DO
790 END DO
791# endif
792# endif
793
794# undef IR_RANGE
795# undef IU_RANGE
796# undef JR_RANGE
797# undef JV_RANGE
798
799 RETURN
subroutine mp_gather2d(ng, model, lbi, ubi, lbj, ubj, tindex, gtype, ascl, amask, a, npts, awrk, setfillval)
integer stdout
integer, dimension(:), allocatable nxyp
integer, dimension(:), allocatable riub
integer, dimension(:), allocatable rjlb
integer, dimension(:), allocatable nxyu
integer, dimension(:), allocatable vjlb
integer isvvel
integer isvbar
integer, dimension(:), allocatable uilb
integer, dimension(:), allocatable nwaterv
integer isvstr
integer, dimension(:), allocatable nxyv
integer, dimension(:), allocatable vilb
integer, dimension(:), allocatable nwateru
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer, dimension(:), allocatable vjub
integer, dimension(:), allocatable rjub
integer, dimension(:), allocatable ujlb
integer isustr
integer, dimension(:), allocatable ujub
integer isubar
integer, dimension(:), allocatable istsur
integer, dimension(:), allocatable nwaterr
integer, dimension(:), allocatable rilb
integer, dimension(:), allocatable uiub
integer, dimension(:), allocatable viub
integer, dimension(:), allocatable nxyr
integer block_count
integer numthreads
logical master
integer, dimension(:), allocatable nstate
Definition mod_param.F:645
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable mstate
Definition mod_param.F:644
integer, dimension(:), allocatable nstr
Definition mod_param.F:646
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
type(t_iobounds), dimension(:), allocatable iobounds
Definition mod_param.F:282
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 ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
integer, dimension(:), allocatable nend
Definition mod_param.F:647
integer, parameter p2dvar
Definition mod_param.F:716
integer, dimension(:), allocatable nt
Definition mod_param.F:489
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
type(io_desc_t), dimension(:), pointer iodesc_dp_resid
type(io_desc_t), dimension(:), pointer iodesc_sp_bvec
type(iosystem_desc_t), dimension(:,:), allocatable, target piosystem
type(io_desc_t), dimension(:), pointer iodesc_dp_bvec
type(io_desc_t), dimension(:), pointer iodesc_sp_sworkd
type(io_desc_t), dimension(:), pointer iodesc_sp_resid
type(io_desc_t), dimension(:), pointer iodesc_dp_sworkd
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
type(t_scalars), dimension(:), allocatable scalars
Definition mod_scalars.F:65
subroutine, public state_iodecomp(ng, iosystem, iotype, iodesc, iovname, ndims)
Definition set_pio.F:1584

References mod_parallel::block_count, mod_param::domain, mod_scalars::ewperiodic, mod_param::iobounds, mod_pio_netcdf::iodesc_dp_bvec, mod_pio_netcdf::iodesc_dp_resid, mod_pio_netcdf::iodesc_dp_sworkd, mod_pio_netcdf::iodesc_sp_bvec, mod_pio_netcdf::iodesc_sp_resid, mod_pio_netcdf::iodesc_sp_sworkd, mod_pio_netcdf::ipioroms, mod_ncparam::isfsur, mod_ncparam::istsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isustr, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvstr, mod_ncparam::isvvel, mod_param::lm, mod_parallel::master, mod_param::mm, distribute_mod::mp_gather2d(), mod_param::mstate, mod_parallel::myrank, mod_param::n, mod_param::nend, mod_param::ngrids, mod_scalars::nsperiodic, mod_param::nstate, mod_param::nstr, mod_param::nt, mod_param::ntilee, mod_param::ntilex, mod_parallel::numthreads, mod_ncparam::nwaterr, mod_ncparam::nwateru, mod_ncparam::nwaterv, mod_ncparam::nxyp, mod_ncparam::nxyr, mod_ncparam::nxyu, mod_ncparam::nxyv, mod_param::p2dvar, mod_pio_netcdf::piosystem, mod_param::r2dvar, mod_ncparam::rilb, mod_ncparam::riub, mod_ncparam::rjlb, mod_ncparam::rjub, mod_scalars::scalars, set_pio_mod::state_iodecomp(), mod_iounits::stdout, mod_param::u2dvar, mod_ncparam::uilb, mod_ncparam::uiub, mod_ncparam::ujlb, mod_ncparam::ujub, mod_param::v2dvar, mod_ncparam::vilb, mod_ncparam::viub, mod_ncparam::vjlb, and mod_ncparam::vjub.

Referenced by wpoints().

Here is the call graph for this function:
Here is the caller graph for this function: