#include "cppdefs.h" #ifdef FULL_GRID # define IR_RANGE IstrT,IendT # define IU_RANGE IstrP,IendT # define JR_RANGE JstrT,JendT # define JV_RANGE JstrP,JendT #else # define IR_RANGE Istr,Iend # define IU_RANGE IstrU,Iend # define JR_RANGE Jstr,Jend # define JV_RANGE JstrV,Jend #endif MODULE wpoints_mod #if defined PROPAGATOR || \ (defined MASKING && (defined READ_WATER || defined WRITE_WATER)) ! !git $Id$ !svn $Id: wpoints.F 1151 2023-02-09 03:08:53Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! This routine sets variables associated with reading and writing of ! ! water points data and/or determining size of packed state vector. ! ! ! !======================================================================= ! implicit none PRIVATE PUBLIC :: wpoints CONTAINS ! !*********************************************************************** SUBROUTINE wpoints (ng, tile, model) !*********************************************************************** ! USE mod_param USE mod_grid ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, model ! ! Local variable declarations. ! # include "tile.h" ! CALL wpoints_tile (ng, tile, model, & & LBi, UBi, LBj, UBj, & & IminS, ImaxS, JminS, JmaxS, & # ifdef MASKING & GRID(ng) % pmask_full, & & GRID(ng) % rmask_full, & & GRID(ng) % umask_full, & & GRID(ng) % vmask_full, & # ifdef PROPAGATOR & GRID(ng) % IJwaterR, & & GRID(ng) % IJwaterU, & & GRID(ng) % IJwaterV, & # endif # endif & GRID(ng) % h) RETURN END SUBROUTINE wpoints ! !*********************************************************************** SUBROUTINE wpoints_tile (ng, tile, model, & & LBi, UBi, LBj, UBj, & & IminS, ImaxS, JminS, JmaxS, & # ifdef MASKING & pmask_full, rmask_full, & & umask_full, vmask_full, & # ifdef PROPAGATOR & IJwaterR, IJwaterU, IJwaterV, & # endif # endif & h) !*********************************************************************** ! USE mod_param USE mod_parallel USE mod_iounits USE mod_ncparam # if defined PROPAGATOR && defined CHECKPOINTING && \ defined PIO_LIB && defined DISTRIBUTE USE mod_pio_netcdf # endif USE mod_scalars # ifdef PROPAGATOR USE mod_storage # endif # ifdef DISTRIBUTE ! USE distribute_mod, ONLY : mp_bcastf, mp_gather2d, mp_reduce # if defined PROPAGATOR && defined CHECKPOINTING && defined PIO_LIB USE set_pio_mod, ONLY : state_iodecomp # endif # endif ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, model integer, intent(in) :: LBi, UBi, LBj, UBj integer, intent(in) :: IminS, ImaxS, JminS, JmaxS ! # ifdef ASSUMED_SHAPE # ifdef MASKING # ifdef PROPAGATOR integer, intent(out) :: IJwaterR(LBi:,LBj:) integer, intent(out) :: IJwaterU(LBi:,LBj:) integer, intent(out) :: IJwaterV(LBi:,LBj:) # endif real(r8), intent(in) :: rmask_full(LBi:,LBj:) real(r8), intent(in) :: pmask_full(LBi:,LBj:) real(r8), intent(in) :: umask_full(LBi:,LBj:) real(r8), intent(in) :: vmask_full(LBi:,LBj:) # endif real(r8), intent(in) :: h(LBi:,LBj:) # else # ifdef MASKING # ifdef PROPAGATOR integer, intent(out) :: IJwaterR(LBi:UBi,LBj:UBj) integer, intent(out) :: IJwaterU(LBi:UBi,LBj:UBj) integer, intent(out) :: IJwaterV(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: rmask_full(LBi:UBi,LBj:UBj) real(r8), intent(in) :: pmask_full(LBi:UBi,LBj:UBj) real(r8), intent(in) :: umask_full(LBi:UBi,LBj:UBj) real(r8), intent(in) :: vmask_full(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: h(LBi:UBi,LBj:UBj) # endif ! ! Local variable declarations. ! # if defined PROPAGATOR && defined CHECKPOINTING && \ defined PIO_LIB && defined DISTRIBUTE logical, save :: first = .TRUE. ! # endif integer :: my_Nxyp, my_Nxyr, my_Nxyu, my_Nxyv # ifdef PROPAGATOR integer :: my_NwaterR, my_NwaterU, my_NwaterV # endif # ifdef PROPAGATOR integer :: Imin, Imax, Jmin, Jmax # endif integer :: Uoff, Voff integer :: NSUB, Npts, i, ic, ij, j # ifdef DISTRIBUTE # ifdef PROPAGATOR integer :: block_size real(r8), dimension(3) :: wp_buffer character (len=3), dimension(3) :: wp_handle # endif # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) real(r8), dimension(4) :: io_buffer character (len=3), dimension(4) :: io_handle # endif # endif # if defined READ_WATER || defined PROPAGATOR real(r8), dimension((Lm(ng)+2)*(Mm(ng)+2)) :: mask # endif # include "set_bounds.h" # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) ! !----------------------------------------------------------------------- ! Determine number of water points. !----------------------------------------------------------------------- ! ! Determine interior points for IO purposes. ! my_Nxyr=0 my_Nxyp=0 my_Nxyu=0 my_Nxyv=0 DO j=Jstr,Jend DO i=Istr,Iend IF (pmask_full(i,j).gt.0.0_r8) my_Nxyp=my_Nxyp+1 IF (rmask_full(i,j).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(i,j).gt.0.0_r8) my_Nxyu=my_Nxyu+1 IF (vmask_full(i,j).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END DO END DO ! ! Determine boundary points for IO purposes. The assigment ! below is done to account for periodicity. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN DO j=Jstr,Jend IF (rmask_full(Istr-1,j).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (vmask_full(Istr-1,j).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END DO END IF IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN DO j=Jstr,Jend IF (pmask_full(Iend+1,j).gt.0.0_r8) my_Nxyp=my_Nxyp+1 IF (rmask_full(Iend+1,j).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(Iend+1,j).gt.0.0_r8) my_Nxyu=my_Nxyu+1 IF (vmask_full(Iend+1,j).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END DO END IF IF (DOMAIN(ng)%Southern_Edge(tile)) THEN DO i=Istr,Iend IF (rmask_full(i,Jstr-1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(i,Jstr-1).gt.0.0_r8) my_Nxyu=my_Nxyu+1 END DO END IF IF (DOMAIN(ng)%Northern_Edge(tile)) THEN DO i=Istr,Iend IF (pmask_full(i,Jend+1).gt.0.0_r8) my_Nxyp=my_Nxyp+1 IF (rmask_full(i,Jend+1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(i,Jend+1).gt.0.0_r8) my_Nxyu=my_Nxyu+1 IF (vmask_full(i,Jend+1).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END DO END IF IF (DOMAIN(ng)%SouthWest_Corner(tile)) THEN IF (rmask_full(Istr-1,Jstr-1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 END IF IF (DOMAIN(ng)%SouthEast_Corner(tile)) THEN IF (rmask_full(Iend+1,Jstr-1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(Iend+1,Jstr-1).gt.0.0_r8) my_Nxyu=my_Nxyu+1 END IF IF (DOMAIN(ng)%NorthWest_Corner(tile)) THEN IF (rmask_full(Istr-1,Jend+1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (vmask_full(Istr-1,Jend+1).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END IF IF (DOMAIN(ng)%NorthEast_Corner(tile)) THEN IF (pmask_full(Iend+1,Jend+1).gt.0.0_r8) my_Nxyp=my_Nxyp+1 IF (rmask_full(Iend+1,Jend+1).gt.0.0_r8) my_Nxyr=my_Nxyr+1 IF (umask_full(Iend+1,Jend+1).gt.0.0_r8) my_Nxyu=my_Nxyu+1 IF (vmask_full(Iend+1,Jend+1).gt.0.0_r8) my_Nxyv=my_Nxyv+1 END IF # endif # ifdef PROPAGATOR ! !----------------------------------------------------------------------- ! Determine number of water points for propagator state vector. !----------------------------------------------------------------------- ! my_NwaterR=0 my_NwaterU=0 my_NwaterV=0 DO j=JR_RANGE DO i=IR_RANGE # ifdef MASKING IF (rmask_full(i,j).gt.0.0_r8) THEN my_NwaterR=my_NwaterR+1 END IF # else my_NwaterR=my_NwaterR+1 # endif END DO DO i=IU_RANGE # ifdef MASKING IF (umask_full(i,j).gt.0.0_r8) THEN my_NwaterU=my_NwaterU+1 END IF # else my_NwaterU=my_NwaterU+1 # endif END DO END DO DO j=JV_RANGE DO i=IR_RANGE # ifdef MASKING IF (vmask_full(i,j).gt.0.0_r8) THEN my_NwaterV=my_NwaterV+1 END IF # else my_NwaterV=my_NwaterV+1 # endif END DO END DO # endif ! ! Determine global number of points (parallel reduction). ! # ifdef DISTRIBUTE NSUB=1 ! distributed-memory # else IF (DOMAIN(ng)%SouthWest_Corner(tile).and. & & DOMAIN(ng)%NorthEast_Corner(tile)) THEN NSUB=1 ! non-tiled application ELSE NSUB=NtileX(ng)*NtileE(ng) ! tiled application END IF # endif !$OMP CRITICAL (MAX_WATER) IF (block_count.eq.0) THEN # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) Nxyp(ng)=0 Nxyr(ng)=0 Nxyu(ng)=0 Nxyv(ng)=0 # endif # ifdef PROPAGATOR NwaterR(ng)=0 NwaterU(ng)=0 NwaterV(ng)=0 # endif END IF # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) Nxyp(ng)=Nxyp(ng)+my_Nxyp Nxyr(ng)=Nxyr(ng)+my_Nxyr Nxyu(ng)=Nxyu(ng)+my_Nxyu Nxyv(ng)=Nxyv(ng)+my_Nxyv # endif # ifdef PROPAGATOR NwaterR(ng)=NwaterR(ng)+my_NwaterR NwaterU(ng)=NwaterU(ng)+my_NwaterU NwaterV(ng)=NwaterV(ng)+my_NwaterV # endif block_count=block_count+1 IF (block_count.eq.NSUB) THEN block_count=0 # ifdef DISTRIBUTE # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) io_buffer(1)=REAL(Nxyp(ng),r8) io_buffer(2)=REAL(Nxyr(ng),r8) io_buffer(3)=REAL(Nxyu(ng),r8) io_buffer(4)=REAL(Nxyv(ng),r8) io_handle(1)='SUM' io_handle(2)='SUM' io_handle(3)='SUM' io_handle(4)='SUM' CALL mp_reduce (ng, model, 4, io_buffer, io_handle) Nxyp(ng)=INT(io_buffer(1)) Nxyr(ng)=INT(io_buffer(2)) Nxyu(ng)=INT(io_buffer(3)) Nxyv(ng)=INT(io_buffer(4)) # endif # ifdef PROPAGATOR wp_buffer(1)=REAL(NwaterR(ng),r8) wp_buffer(2)=REAL(NwaterU(ng),r8) wp_buffer(3)=REAL(NwaterV(ng),r8) wp_handle(1)='SUM' wp_handle(2)='SUM' wp_handle(3)='SUM' CALL mp_reduce (ng, model, 3, wp_buffer, wp_handle) NwaterR(ng)=INT(wp_buffer(1)) NwaterU(ng)=INT(wp_buffer(2)) NwaterV(ng)=INT(wp_buffer(3)) # endif # endif # if defined MASKING && (defined READ_WATER || defined WRITE_WATER) IOBOUNDS(ng) % xy_psi = Nxyp(ng) IOBOUNDS(ng) % xy_rho = Nxyr(ng) IOBOUNDS(ng) % xy_u = Nxyu(ng) IOBOUNDS(ng) % xy_v = Nxyv(ng) # endif END IF !$OMP END CRITICAL (MAX_WATER) # ifdef PROPAGATOR ! !----------------------------------------------------------------------- ! Determine size of state vector. !----------------------------------------------------------------------- # if defined FORCING_SV || \ defined SO_SEMI || defined STOCHASTIC_OPT ! ! The state vector may contain free-surface, momentum, tracers, ! surface momentum stresses, amd surface tracer fluxes. ! Mstate(ng)=0 IF (SCALARS(ng)%Fstate(isFsur)) THEN Mstate(ng)=Mstate(ng)+NwaterR(ng) END IF # ifdef SOLVE3D IF (SCALARS(ng)%Fstate(isUvel)) THEN Mstate(ng)=Mstate(ng)+N(ng)*NwaterU(ng) END IF IF (SCALARS(ng)%Fstate(isVvel)) THEN Mstate(ng)=Mstate(ng)+N(ng)*NwaterV(ng) END IF DO i=1,NT(ng) IF (SCALARS(ng)%Fstate(isTvar(i))) THEN Mstate(ng)=Mstate(ng)+N(ng)*NwaterR(ng) END IF END DO # else IF (SCALARS(ng)%Fstate(isUbar)) THEN Mstate(ng)=Mstate(ng)+NwaterU(ng) END IF IF (SCALARS(ng)%Fstate(isVbar)) THEN Mstate(ng)=Mstate(ng)+NwaterV(ng) END IF # endif IF (SCALARS(ng)%Fstate(isUstr)) THEN Mstate(ng)=Mstate(ng)+NwaterU(ng) END IF IF (SCALARS(ng)%Fstate(isVstr)) THEN Mstate(ng)=Mstate(ng)+NwaterV(ng) END IF # ifdef SOLVE3D DO i=1,NT(ng) IF (SCALARS(ng)%Fstate(isTsur(i))) THEN Mstate(ng)=Mstate(ng)+NwaterR(ng) END IF END DO # endif # else # ifdef SOLVE3D ! ! The state vector contains free-surface, 3D momentum components, and ! tracers. ! Mstate(ng)=NwaterR(ng)+ & & NwaterU(ng)*N(ng)+ & & NwaterV(ng)*N(ng)+ & & NwaterR(ng)*N(ng)*NT(ng) # else ! ! The state vector contains free-surface and 2D momentum components. ! Mstate(ng)=NwaterR(ng)+ & & NwaterU(ng)+ & & NwaterV(ng) # endif # endif # ifdef DISTRIBUTE ! ! Deternine state vector starting and ending blocking indices. ! Nstate(ng)=(Mstate(ng)+numthreads-1)/numthreads block_size=Nstate(ng) IF (Master) THEN WRITE (stdout,10) ng, Mstate(ng), Nstate(ng) 10 FORMAT (/,' Size of FULL state arrays for propagator of grid ', & & i2.2,', Mstate = ', i10,/,38x, & & 'NODE partition, Nstate = ', i10,/,/, & & 9x,'Node',7x,'Nstr',7x,'Nend',7x,'Size',/) DO ic=0,numthreads-1 i=1+ic*block_size j=MIN(Mstate(ng),i+block_size-1) WRITE (stdout,20) ic, i, j, j-i+1 20 FORMAT (11x, i2, 3(1x,i10)) END DO WRITE (stdout,'(/)') END IF Nstr(ng)=1+MyRank*block_size Nend(ng)=Nstr(ng)+block_size-1 Mstate(ng)=Nstate(ng)*numthreads # else Nstr(ng)=1 Nend(ng)=Mstate(ng) Nstate(ng)=Mstate(ng) IF (Master) THEN WRITE (stdout,10) ng, Mstate(ng) 10 FORMAT (' Size of FULL state arrays for propagator of grid ', & & i2.2, ', Mstate = ', i10,/) END IF # endif # endif # if defined PIO_LIB && defined DISTRIBUTE && defined CHECKPOINTING ! !----------------------------------------------------------------------- ! Set I/O decomposition descriptors for single and double precision ! state propagator data !----------------------------------------------------------------------- ! IF (first) THEN IF (ng.eq.Ngrids) first=.FALSE. ! CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_real, & & ioDesc_sp_Bvec(ng), & & 'Bvec', 2) CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_real, & & ioDesc_sp_resid(ng), & & 'resid', 1) CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_real, & & ioDesc_sp_SworkD(ng), & & 'SworkD', 1) ! CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_double, & & ioDesc_dp_Bvec(ng), & & 'Bvec', 2) CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_double, & & ioDesc_dp_resid(ng), & & 'resid', 1) CALL state_iodecomp (ng, pioSystem(IpioROMS,ng), PIO_double, & & ioDesc_dp_SworkD(ng), & & 'SworkD', 1) END IF # endif # if ((defined READ_WATER && defined DISTRIBUTE) || \ defined PROPAGATOR) && defined MASKING ! !----------------------------------------------------------------------- ! If distributed-memory, set process water indices arrays. !----------------------------------------------------------------------- ! ! Set offsets for U- and V-variables due to periodic boundary ! conditions. Recall that in East-West periodic boundary conditions ! IstrU=1. Otherwise, IstrU=2. Similarly, in North-South periodic ! applications IstrV=1 or else IstrV=2. ! IF (EWperiodic(ng)) THEN Uoff=0 ELSE Uoff=1 END IF ! IF (NSperiodic(ng)) THEN Voff=0 ELSE Voff=1 END IF # ifdef DISTRIBUTE ! ! Gather masking at RHO-points from all nodes. ! CALL mp_gather2d (ng, model, LBi, UBi, LBj, UBj, & & 0, r2dvar, 1.0_r8, & & rmask_full, rmask_full, Npts, mask, .FALSE.) CALL mp_bcastf (ng, model, mask) # else ! ! Load masking at RHO-points in mask vector. ! ic=0 DO j=0,Mm(ng)+1 DO i=0,Lm(ng)+1 ic=ic+1 mask(ic)=rmask_full(i,j) END DO END DO Npts=ic # endif # if defined READ_WATER && defined DISTRIBUTE ! ! Determine water RHO-points for IO purposes. ! ic=0 DO ij=1,Npts IF (mask(ij).gt.0.0_r8) THEN ic=ic+1 SCALARS(ng)%IJwater(ic,r2dvar)=ij END IF END DO my_Nxyr=ic # endif # ifdef PROPAGATOR ! ! Determine water points for the propagator state vector at RHO-points. ! ic=0 ij=0 # ifdef FULL_GRID Imin=0 Imax=Lm(ng)+1 Jmin=0 Jmax=Mm(ng)+1 # else Imin=1 Imax=Lm(ng) Jmin=1 Jmax=Mm(ng) # endif DO j=0,Mm(ng)+1 DO i=0,Lm(ng)+1 ij=ij+1 IF ((mask(ij).gt.0.0_r8).and. & & (Imin.le.i).and.(i.le.Imax).and. & & (Jmin.le.j).and.(j.le.Jmax)) THEN ic=ic+1 mask(ij)=REAL(ic,r8) ELSE mask(ij)=0.0_r8 END IF END DO END DO ! ij=0 DO j=0,Mm(ng)+1 DO i=0,Lm(ng)+1 ij=ij+1 IF ((rILB(ng).le.i).and.(i.le.rIUB(ng)).and. & & (rJLB(ng).le.j).and.(j.le.rJUB(ng))) THEN IF (mask(ij).gt.0.0_r8) THEN IJwaterR(i,j)=INT(mask(ij)) ELSE IJwaterR(i,j)=0 END IF END IF END DO END DO # endif # if defined READ_WATER && defined DISTRIBUTE ! ! Gather masking at PSI-points from all nodes. ! CALL mp_gather2d (ng, model, LBi, UBi, LBj, UBj, & & 0, p2dvar, 1.0_r8, & & pmask_full, pmask_full, Npts, mask, .FALSE.) CALL mp_bcastf (ng, model, mask) ! ! Determine water PSI-points for IO purposes. ! ic=0 DO ij=1,Npts IF (mask(ij).gt.0.0_r8) THEN ic=ic+1 SCALARS(ng)%IJwater(ic,p2dvar)=ij END IF END DO my_Nxyp=ic # endif # ifdef DISTRIBUTE ! ! Gather masking at U-points from all nodes. ! CALL mp_gather2d (ng, model, LBi, UBi, LBj, UBj, & & 0, u2dvar, 1.0_r8, & & umask_full, umask_full, Npts, mask, .FALSE.) CALL mp_bcastf (ng, model, mask) # else ! ! Load masking at U-points in mask vector. ! ic=0 DO j=0,Mm(ng)+1 DO i=1,Lm(ng)+1 ic=ic+1 mask(ic)=umask_full(i,j) END DO END DO Npts=ic # endif # if defined READ_WATER && defined DISTRIBUTE ! ! Determine water U-points for IO purposes. ! ic=0 DO ij=1,Npts IF (mask(ij).gt.0.0_r8) THEN ic=ic+1 SCALARS(ng)%IJwater(ic,u2dvar)=ij END IF END DO my_Nxyu=ic # endif # ifdef PROPAGATOR ! ! Determine water points for the propagator state vector at U-points. ! ic=0 ij=0 # ifdef FULL_GRID Imin=1 Imax=Lm(ng)+1 Jmin=0 Jmax=Mm(ng)+1 # else Imin=1+Uoff Imax=Lm(ng) Jmin=1 Jmax=Mm(ng) # endif DO j=0,Mm(ng)+1 DO i=1,Lm(ng)+1 ij=ij+1 IF ((mask(ij).gt.0.0_r8).and. & & (Imin.le.i).and.(i.le.Imax).and. & & (Jmin.le.j).and.(j.le.Jmax)) THEN ic=ic+1 mask(ij)=REAL(ic,r8) ELSE mask(ij)=0.0_r8 END IF END DO END DO ! ij=0 DO j=0,Mm(ng)+1 DO i=1,Lm(ng)+1 ij=ij+1 IF ((uILB(ng).le.i).and.(i.le.uIUB(ng)).and. & & (uJLB(ng).le.j).and.(j.le.uJUB(ng))) THEN IF (mask(ij).gt.0.0_r8) THEN IJwaterU(i,j)=INT(mask(ij)) ELSE IJwaterU(i,j)=0 END IF END IF END DO END DO # endif # ifdef DISTRIBUTE ! ! Gather masking at V-points from all nodes. ! CALL mp_gather2d (ng, model, LBi, UBi, LBj, UBj, & & 0, v2dvar, 1.0_r8, & & vmask_full, vmask_full, Npts, mask, .FALSE.) CALL mp_bcastf (ng, model, mask) # else ! ! Load masking at V-points in mask vector. ! ic=0 DO j=1,Mm(ng)+1 DO i=0,Lm(ng)+1 ic=ic+1 mask(ic)=vmask_full(i,j) END DO END DO Npts=ic # endif # if defined READ_WATER && defined DISTRIBUTE ! ! Determine water V-points for IO purposes. ! ic=0 DO ij=1,Npts IF (mask(ij).gt.0.0_r8) THEN ic=ic+1 SCALARS(ng)%IJwater(ic,v2dvar)=ij END IF END DO my_Nxyv=ic # endif # ifdef PROPAGATOR ! ! Determine water points for the propagator state vector at V-points. ! ic=0 ij=0 # ifdef FULL_GRID Imin=0 Imax=Lm(ng)+1 Jmin=1 Jmax=Mm(ng)+1 # else Imin=1 Imax=Lm(ng) Jmin=1+Voff Jmax=Mm(ng) # endif DO j=1,Mm(ng)+1 DO i=0,Lm(ng)+1 ij=ij+1 IF ((mask(ij).gt.0.0_r8).and. & & (Imin.le.i).and.(i.le.Imax).and. & & (Jmin.le.j).and.(j.le.Jmax)) THEN ic=ic+1 mask(ij)=REAL(ic,r8) ELSE mask(ij)=0.0_r8 END IF END DO END DO ! ij=0 DO j=1,Mm(ng)+1 DO i=0,Lm(ng)+1 ij=ij+1 IF ((vILB(ng).le.i).and.(i.le.vIUB(ng)).and. & & (vJLB(ng).le.j).and.(j.le.vJUB(ng))) THEN IF (mask(ij).gt.0.0_r8) THEN IJwaterV(i,j)=INT(mask(ij)) ELSE IJwaterV(i,j)=0 END IF END IF END DO END DO # endif # endif # undef IR_RANGE # undef IU_RANGE # undef JR_RANGE # undef JV_RANGE RETURN END SUBROUTINE wpoints_tile #endif END MODULE wpoints_mod