SUBROUTINE get_bounds (ng, tile, gtype, Nghost, Itile, Jtile, & & LBi, UBi, LBj, UBj) ! !git $Id$ !svn $Id: get_bounds.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 compute grid bounds in the I- and J-directions. ! ! ! ! On Input: ! ! ! ! ng Nested grid number. ! ! tile Domain partition. ! ! gtype C-grid type. If zero, compute array allocation bounds.! ! Otherwise, compute bounds for IO processing. ! ! Nghost Number of ghost-points in the halo region: ! ! Nghost = 0, compute non-overlapping bounds. ! ! Nghost > 0, compute overlapping bounds. ! ! ! ! On Output: ! ! ! ! Itile I-tile coordinate (a value from 0 to NtileI(ng)). ! ! Jtile J-tile coordinate (a value from 0 to NtileJ(ng)). ! ! LBi I-dimension Lower bound. ! ! UBi I-dimension Upper bound. ! ! LBj J-dimension Lower bound. ! ! UBj J-dimension Upper bound. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_ncparam USE mod_scalars ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, gtype, Nghost integer, intent(out) :: Itile, Jtile, LBi, UBi, LBj, UBj ! ! Local variable declarations. ! integer :: Imin, Imax, Jmin, Jmax integer :: Iend, Istr, Jend, Jstr integer :: IstrM, IstrR, IstrU, IendR integer :: JstrM, JstrR, JstrV, JendR integer :: IstrB, IendB, IstrP, IendP, IstrT, IendT integer :: JstrB, JendB, JstrP, JendP, JstrT, JendT integer :: Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1 integer :: Iendp1, Iendp2, Iendp2i, Iendp3 integer :: Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1 integer :: Jendp1, Jendp2, Jendp2i, Jendp3 integer :: MyType ! !----------------------------------------------------------------------- ! Set array bounds in the I- and J-direction for distributed-memory ! configurations. !----------------------------------------------------------------------- ! ! Set first and last grid-points according to staggered C-grid ! classification. If gtype = 0, it returns the values needed for ! array allocation. Otherwise, it returns the values needed for IO ! processing. ! MyType=ABS(gtype) IF (MyType.eq.0) THEN IF (EWperiodic(ng)) THEN IF (NSperiodic(ng)) THEN Imin=-NghostPoints Imax=Im(ng)+NghostPoints Jmin=-NghostPoints Jmax=Jm(ng)+NghostPoints ELSE Imin=-NghostPoints Imax=Im(ng)+NghostPoints Jmin=0 Jmax=Jm(ng)+1 END IF ELSE IF (NSperiodic(ng)) THEN Imin=0 Imax=Im(ng)+1 Jmin=-NghostPoints Jmax=Jm(ng)+NghostPoints ELSE Imin=0 Imax=Im(ng)+1 Jmin=0 Jmax=Jm(ng)+1 END IF END IF ELSE IF ((MyType.eq.p2dvar).or.(MyType.eq.u2dvar).or. & & (MyType.eq.p3dvar).or.(MyType.eq.u3dvar)) THEN Imin=1 ELSE Imin=0 END IF Imax=Lm(ng)+1 IF ((MyType.eq.p2dvar).or.(MyType.eq.v2dvar).or. & & (MyType.eq.p3dvar).or.(MyType.eq.v3dvar)) THEN Jmin=1 ELSE Jmin=0 END IF Jmax=Mm(ng)+1 END IF ! ! Set physical, overlapping (Nghost>0) or non-overlapping (Nghost=0) ! grid bounds according to tile rank. ! CALL get_tile (ng, tile, & & Itile, Jtile, & & Istr, Iend, Jstr, Jend, & & IstrM, IstrR, IstrU, IendR, & & JstrM, JstrR, JstrV, JendR, & & IstrB, IendB, IstrP, IendP, IstrT, IendT, & & JstrB, JendB, JstrP, JendP, JstrT, JendT, & & Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1, & & Iendp1, Iendp2, Iendp2i, Iendp3, & & Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1, & & Jendp1, Jendp2, Jendp2i, Jendp3) ! IF ((Itile.eq.-1).or.(Itile.eq.0)) THEN LBi=Imin ELSE LBi=Istr-Nghost END IF IF ((Itile.eq.-1).or.(Itile.eq.(NtileI(ng)-1))) THEN UBi=Imax ELSE UBi=Iend+Nghost END IF IF ((Jtile.eq.-1).or.(Jtile.eq.0)) THEN LBj=Jmin ELSE LBj=Jstr-Nghost END IF IF ((Jtile.eq.-1).or.(Jtile.eq.(NtileJ(ng)-1))) THEN UBj=Jmax ELSE UBj=Jend+Nghost END IF RETURN END SUBROUTINE get_bounds SUBROUTINE get_domain (ng, tile, gtype, Nghost, & & epsilon, Lfullgrid, & & Xmin, Xmax, Ymin, Ymax) ! !======================================================================= ! ! ! This routine computes tile minimum and maximum fractional grid ! ! coordinates. ! ! ! ! On Input: ! ! ! ! ng Nested grid number. ! ! tile Domain partition. ! ! Nghost Number of ghost-points in the halo region: ! ! Nghost = 0, compute non-overlapping coordinates. ! ! Nghost > 0, compute overlapping bounds. ! ! gtype C-grid type ! ! epsilon Small value to add to Xmax and Ymax when the tile ! ! is lying on the eastern and northern boundaries ! ! of the grid. This is usefull when processing ! ! observations. ! ! Lfullgrid Switch to include interior and boundaries points ! ! (TRUE) or just interior points (FALSE). ! ! ! ! On Output: ! ! ! ! Xmin Minimum tile fractional X-coordinate. ! ! Xmax Maximum tile fractional X-coordinate. ! ! Ymin Minimum tile fractional Y-coordinate. ! ! Ymax Maximum tile fractional Y-coordinate. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_ncparam ! implicit none ! ! Imported variable declarations. ! logical, intent(in) :: Lfullgrid integer, intent(in) :: ng, tile, gtype, Nghost real(r8), intent(in) :: epsilon real(r8), intent(out) :: Xmin, Xmax, Ymin, Ymax ! ! Local variable declarations. ! integer :: Imin, Imax, Jmin, Jmax integer :: Itile, Jtile ! !----------------------------------------------------------------------- ! Computes tile minimum and maximum fractional-grid coordinates. !----------------------------------------------------------------------- ! CALL get_bounds (ng, tile, gtype, Nghost, Itile, Jtile, & & Imin, Imax, Jmin, Jmax) ! ! Include interior and boundary points. ! IF (Lfullgrid) THEN IF ((Itile.eq.0).and. & & ((gtype.eq.r2dvar).or.(gtype.eq.r3dvar).or. & & (gtype.eq.v2dvar).or.(gtype.eq.v3dvar))) THEN Xmin=REAL(Imin,r8) ELSE Xmin=REAL(Imin,r8)-0.5_r8 END IF IF (Itile.eq.(NtileI(ng)-1)) THEN IF ((gtype.eq.u2dvar).or.(gtype.eq.u3dvar)) THEN Xmax=REAL(Imax,r8)-0.5_r8 ELSE Xmax=REAL(Imax,r8) END IF ELSE Xmax=REAL(Imax,r8)+0.5_r8 END IF IF ((Jtile.eq.0).and. & & ((gtype.eq.r2dvar).or.(gtype.eq.r3dvar).or. & & (gtype.eq.u2dvar).or.(gtype.eq.u3dvar))) THEN Ymin=REAL(Jmin,r8) ELSE Ymin=REAL(Jmin,r8)-0.5_r8 END IF IF (Jtile.eq.(NtileJ(ng)-1)) THEN IF ((gtype.eq.v2dvar).or.(gtype.eq.v3dvar)) THEN Ymax=REAL(Jmax,r8)-0.5_r8 ELSE Ymax=REAL(Jmax,r8) END IF ELSE Ymax=REAL(Jmax,r8)+0.5_r8 END IF ! ! Include only interior points. ! ELSE IF (Itile.eq.0) THEN IF ((gtype.eq.u2dvar).or.(gtype.eq.u3dvar)) THEN Xmin=REAL(Imin,r8) ELSE Xmin=REAL(Imin,r8)+0.5_r8 END IF ELSE Xmin=REAL(Imin,r8)-0.5_r8 END IF IF (Itile.eq.(NtileI(ng)-1)) THEN IF ((gtype.eq.u2dvar).or.(gtype.eq.u3dvar)) THEN Xmax=REAL(Imax,r8)-1.0_r8 ELSE Xmax=REAL(Imax,r8)-0.5_r8 END IF ELSE Xmax=REAL(Imax,r8)+0.5_r8 END IF IF (Jtile.eq.0) THEN IF ((gtype.eq.v2dvar).or.(gtype.eq.v3dvar)) THEN Ymin=REAL(Jmin,r8) ELSE Ymin=REAL(Jmin,r8)+0.5 END IF ELSE Ymin=REAL(Jmin,r8)-0.5_r8 END IF IF (Jtile.eq.(NtileJ(ng)-1)) THEN IF ((gtype.eq.v2dvar).or.(gtype.eq.v3dvar)) THEN Ymax=REAL(Jmax,r8)-1.0_r8 ELSE Ymax=REAL(Jmax,r8)-0.5_r8 END IF ELSE Ymax=REAL(Jmax,r8)+0.5_r8 END IF END IF ! ! If tile lie at the grid eastern or northen boundary, add provided ! offset value to allow processing at those boundaries. ! IF (Itile.eq.(NtileI(ng)-1)) THEN Xmax=Xmax+epsilon END IF IF (Jtile.eq.(NtileJ(ng)-1)) THEN Ymax=Ymax+epsilon END IF RETURN END SUBROUTINE get_domain SUBROUTINE get_domain_edges (ng, tile, & & Eastern_Edge, & & Western_Edge, & & Northern_Edge, & & Southern_Edge, & & NorthEast_Corner, & & NorthWest_Corner, & & SouthEast_Corner, & & SouthWest_Corner, & & NorthEast_Test, & & NorthWest_Test, & & SouthEast_Test, & & SouthWest_Test) ! !======================================================================= ! ! ! This routine sets the logical switches (T/F) needed for processing ! ! model variables in tiles adjacent to the domain boundary edges. It ! ! facilitates complicated nesting configurations. ! ! ! ! On Input: ! ! ! ! ng Nested grid number. ! ! tile Domain partition. ! ! ! ! On Output: ! ! ! ! Eastern_Edge tile next to the domain eastern boundary ! ! Western_Edge tile next to the domain western boundary ! ! Northern_Edge tile next to the domain northern boundary ! ! Southern_Edge tile next to the domain southern boundary ! ! ! ! NorthEast_Corner tile next to the domain northeastern corner ! ! NorthWest_Corner tile next to the domain northwestern corner ! ! SouthEast_Corner tile next to the domain southeastern corner ! ! SouthWest_Corner tile next to the domain southwestern corner ! ! ! ! NorthEast_Test test for tiles in the northeastern corner ! ! NorthWest_Test test for tiles in the northwestern corner ! ! SouthEast_Test test for tiles in the southeastern corner ! ! SouthWest_Test test for tiles in the southwestern corner ! ! ! !======================================================================= ! USE mod_param ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile logical, intent(out) :: Eastern_Edge logical, intent(out) :: Western_Edge logical, intent(out) :: Northern_Edge logical, intent(out) :: Southern_Edge logical, intent(out) :: NorthEast_Corner logical, intent(out) :: NorthWest_Corner logical, intent(out) :: SouthEast_Corner logical, intent(out) :: SouthWest_Corner logical, intent(out) :: NorthEast_Test logical, intent(out) :: NorthWest_Test logical, intent(out) :: SouthEast_Test logical, intent(out) :: SouthWest_Test ! ! Local variable declarations. ! integer :: Istr, Iend, Jstr, Jend integer :: Itile, Jtile ! !----------------------------------------------------------------------- ! Compute Itile and Jtile. !----------------------------------------------------------------------- ! IF (tile.eq.-1) THEN Itile=-1 Jtile=-1 Istr=1 Iend=Lm(ng) Jstr=1 Jend=Mm(ng) ELSE CALL tile_bounds_2d (ng, tile, Lm(ng), Mm(ng), Itile, Jtile, & & Istr, Iend, Jstr, Jend) END IF ! !----------------------------------------------------------------------- ! Set the logical switches (T/F) needed for processing model variables ! in tiles adjacent to the domain boundary edges. !----------------------------------------------------------------------- ! ! HGA: Need to add the logic for composed grids. IF (tile.eq.-1) THEN ! ! Set switches for the full grid (tile=-1) to TRUE, since it contains ! all the boundary edges and corners. This is a special case use for ! other purposes and need only in routine "var_bounds". ! Western_Edge=.TRUE. Eastern_Edge=.TRUE. Southern_Edge=.TRUE. Northern_Edge=.TRUE. SouthWest_Test=.TRUE. SouthEast_Test=.TRUE. NorthWest_Test=.TRUE. NorthEast_Test=.TRUE. SouthWest_Corner=.TRUE. SouthEast_Corner=.TRUE. NorthWest_Corner=.TRUE. NorthEast_Corner=.TRUE. ELSE ! ! Is the tile adjacent to the western or eastern domain ! boundary edge? ! IF (Itile.eq.0) THEN Western_Edge=.TRUE. ! (Istr.eq.1) ELSE Western_Edge=.FALSE. ! elsewhere END IF IF (Itile.eq.(NtileI(ng)-1)) THEN Eastern_Edge=.TRUE. ! (Iend.eq.Lm(ng)) ELSE Eastern_Edge=.FALSE. ! elsewhere END IF ! ! Is the tile adjacent to the southern or northern domain ! boundary edge? ! IF (Jtile.eq.0) THEN Southern_Edge=.TRUE. ! (Jstr.eq.1) ELSE Southern_Edge=.FALSE. ! elsewhere END IF IF (Jtile.eq.(NtileJ(ng)-1)) THEN Northern_Edge=.TRUE. ! (Jend.eq.Mm(ng)) ELSE Northern_Edge=.FALSE. ! elsewhere END IF ! ! Is the tile adjacent to the southwestern domain corner? ! IF ((Itile.eq.0).and. & & (Jtile.eq.0)) THEN SouthWest_Corner=.TRUE. ! (Istr.eq.1).and. SouthWest_Test =.TRUE. ! (Jstr.eq.1) ELSE SouthWest_Corner=.FALSE. ! elsewhere SouthWest_Test =.TRUE. END IF ! ! Is the tile adjacent to the southeastern domain corner? ! IF ((Itile.eq.(NtileI(ng)-1)).and. & & (Jtile.eq.0)) THEN SouthEast_Corner=.TRUE. ! (Iend.eq.Lm(ng)).and. SouthEast_Test =.TRUE. ! (Jstr.eq.1) ELSE SouthEast_Corner=.FALSE. ! elsewhere SouthEast_Test =.TRUE. END IF ! ! Is the tile adjacent to the northwestern domain corner? ! IF ((Itile.eq.0).and. & & (Jtile.eq.(NtileJ(ng)-1))) THEN NorthWest_Corner=.TRUE. ! (Istr.eq.1).and. NorthWest_Test =.TRUE. ! (Jend.eq.Mm(ng)) ELSE NorthWest_Corner=.FALSE. ! elsewhere NorthWest_Test =.TRUE. END IF ! ! Is the tile adjacent to the northeastern domain corner? ! IF ((Itile.eq.(NtileI(ng)-1)).and. & & (Jtile.eq.(NtileJ(ng)-1))) THEN NorthEast_Corner=.TRUE. ! (Iend.eq.Lm(ng)).and. NorthEast_Test =.TRUE. ! (Jend.eq.Mm(ng)) ELSE NorthEast_Corner=.FALSE. ! elsewhere NorthEast_Test =.TRUE. END IF END IF RETURN END SUBROUTINE get_domain_edges SUBROUTINE get_iobounds (ng) ! !======================================================================= ! ! ! This routine computes the horizontal lower bound, upper bound, and ! ! grid size for IO (NetCDF) variables. Nested grids require special ! ! attention due to their connetivity. ! ! ! ! On Input: ! ! ! ! ng Nested grid number. ! ! ! ! On Output, the horizontal lower/upper bounds and grid size for ! ! each variable type and nested grid number are loaded into the ! ! IOBOUNDS structure which is declared in module MOD_PARAM: ! ! ! ! IOBOUNDS(ng) % ILB_psi I-direction lower bound (PSI) ! ! IOBOUNDS(ng) % IUB_psi I-direction upper bound (PSI) ! ! IOBOUNDS(ng) % JLB_psi J-direction lower bound (PSI) ! ! IOBOUNDS(ng) % JUB_psi J-direction upper bound (PSI) ! ! ! ! IOBOUNDS(ng) % ILB_rho I-direction lower bound (RHO) ! ! IOBOUNDS(ng) % IUB_rho I-direction upper bound (RHO) ! ! IOBOUNDS(ng) % JLB_rho J-direction lower bound (RHO) ! ! IOBOUNDS(ng) % JUB_rho J-direction upper bound (RHO) ! ! ! ! IOBOUNDS(ng) % ILB_u I-direction lower bound (U) ! ! IOBOUNDS(ng) % IUB_u I-direction upper bound (U) ! ! IOBOUNDS(ng) % JLB_u J-direction lower bound (U) ! ! IOBOUNDS(ng) % JUB_u J-direction upper bound (U) ! ! ! ! IOBOUNDS(ng) % ILB_v I-direction lower bound (V) ! ! IOBOUNDS(ng) % IUB_v I-direction upper bound (V) ! ! IOBOUNDS(ng) % JLB_v J-direction lower bound (V) ! ! IOBOUNDS(ng) % JUB_v J-direction upper bound (V) ! ! ! ! IOBOUNDS(ng) % xi_psi Number of I-direction points (PSI) ! ! IOBOUNDS(ng) % xi_rho Number of I-direction points (RHO) ! ! IOBOUNDS(ng) % xi_u Number of I-direction points (U) ! ! IOBOUNDS(ng) % xi_v Number of I-direction points (V) ! ! ! ! IOBOUNDS(ng) % eta_psi Number of J-direction points (PSI) ! ! IOBOUNDS(ng) % eta_rho Number of J-direction points (RHO) ! ! IOBOUNDS(ng) % eta_u Number of I-direction points (U) ! ! IOBOUNDS(ng) % eta_v Number of I-direction points (V) ! ! ! !======================================================================= ! USE mod_param ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng ! !----------------------------------------------------------------------- ! Set IO lower/upper bounds and grid size for each C-grid type ! variable. !----------------------------------------------------------------------- ! ! Recall that in non-nested applications the horizontal range, ! including interior and boundary points, for all variable types ! are: ! ! PSI-type [xi_psi, eta_psi] = [1:Lm(ng)+1, 1:Mm(ng)+1] ! RHO-type [xi_rho, eta_rho] = [0:Lm(ng)+1, 0:Mm(ng)+1] ! U-type [xi_u, eta_u ] = [1:Lm(ng)+1, 0:Mm(ng)+1] ! V-type [xi_v, eta_v ] = [0:Lm(ng)+1, 1:Mm(ng)+1] ! IOBOUNDS(ng) % ILB_psi = 1 IOBOUNDS(ng) % IUB_psi = Lm(ng)+1 IOBOUNDS(ng) % JLB_psi = 1 IOBOUNDS(ng) % JUB_psi = Mm(ng)+1 ! IOBOUNDS(ng) % ILB_rho = 0 IOBOUNDS(ng) % IUB_rho = Lm(ng)+1 IOBOUNDS(ng) % JLB_rho = 0 IOBOUNDS(ng) % JUB_rho = Mm(ng)+1 ! IOBOUNDS(ng) % ILB_u = 1 IOBOUNDS(ng) % IUB_u = Lm(ng)+1 IOBOUNDS(ng) % JLB_u = 0 IOBOUNDS(ng) % JUB_u = Mm(ng)+1 ! IOBOUNDS(ng) % ILB_v = 0 IOBOUNDS(ng) % IUB_v = Lm(ng)+1 IOBOUNDS(ng) % JLB_v = 1 IOBOUNDS(ng) % JUB_v = Mm(ng)+1 ! ! Set IO NetCDF files horizontal dimension size. Recall that NetCDF ! does not support arrays with zero index as an array element. ! IOBOUNDS(ng) % IorJ = BOUNDS(ng) % UBij - & & BOUNDS(ng) % LBij + 1 ! IOBOUNDS(ng) % xi_psi = IOBOUNDS(ng) % IUB_psi - & & IOBOUNDS(ng) % ILB_psi + 1 IOBOUNDS(ng) % xi_rho = IOBOUNDS(ng) % IUB_rho - & & IOBOUNDS(ng) % ILB_rho + 1 IOBOUNDS(ng) % xi_u = IOBOUNDS(ng) % IUB_u - & & IOBOUNDS(ng) % ILB_u + 1 IOBOUNDS(ng) % xi_v = IOBOUNDS(ng) % IUB_v - & & IOBOUNDS(ng) % ILB_v + 1 ! IOBOUNDS(ng) % eta_psi = IOBOUNDS(ng) % JUB_psi - & & IOBOUNDS(ng) % JLB_psi + 1 IOBOUNDS(ng) % eta_rho = IOBOUNDS(ng) % JUB_rho - & & IOBOUNDS(ng) % JLB_rho + 1 IOBOUNDS(ng) % eta_u = IOBOUNDS(ng) % JUB_u - & & IOBOUNDS(ng) % JLB_u + 1 IOBOUNDS(ng) % eta_v = IOBOUNDS(ng) % JUB_v - & & IOBOUNDS(ng) % JLB_v + 1 RETURN END SUBROUTINE get_iobounds SUBROUTINE get_tile (ng, tile, & & Itile, Jtile, & & Istr, Iend, Jstr, Jend, & & IstrM, IstrR, IstrU, IendR, & & JstrM, JstrR, JstrV, JendR, & & IstrB, IendB, IstrP, IendP, IstrT, IendT, & & JstrB, JendB, JstrP, JendP, JstrT, JendT, & & Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1, & & Iendp1, Iendp2, Iendp2i, Iendp3, & & Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1, & & Jendp1, Jendp2, Jendp2i, Jendp3) ! !======================================================================= ! ! ! This routine computes the starting and ending horizontal indices ! ! for each sub-domain partition or tile. ! ! ! ! On Input: ! ! ! ! ng Nested grid number (integer). ! ! tile Sub-domain partition. ! ! ! ! On Output: ! ! ! ! Itile I-tile coordinate (a value from 0 to NtileI(ng)). ! ! Jtile J-tile coordinate (a value from 0 to NtileJ(ng)). ! ! ! ! Istr Starting tile index in the I-direction. ! ! Iend Ending tile index in the I-direction. ! ! Jstr Starting tile index in the J-direction. ! ! Jend Ending tile index in the J-direction. ! ! ! ! IstrR Starting tile index in the I-direction (RHO-points) ! ! IstrU Starting tile index in the I-direction (U-points) ! ! IendR Ending tile index in the I-direction (RHO_points) ! ! ! ! JstrR Starting tile index in the J-direction (RHO-points) ! ! JstrV Starting tile index in the J-direction (V-points) ! ! JendR Ending tile index in the J-direction (RHO_points) ! ! ! ! IstrB Starting nest tile in the I-direction (RHO-, V-obc) ! ! IendB Ending nest tile in the I-direction (RHO-, V-obc) ! ! IstrM Starting nest tile in the I-direction (PSI-, U-obc) ! ! IstrP Starting nest tile in the I-direction (PSI, U-points) ! ! IendP Ending nest tile in the I-direction (PSI) ! ! IstrT Starting nest tile in the I-direction (RHO-points) ! ! IendT Ending nest tile in the I-direction (RHO_points) ! ! ! ! JstrB Starting nest tile in the J-direction (RHO-, U-obc) ! ! JendB Ending nest tile in the J-direction (RHO-, U-obc) ! ! JstrM Starting nest tile in the J-direction (PSI-, V-obc) ! ! JstrP Starting nest tile in the J-direction (PSI, V-points) ! ! JendP Ending nest tile in the J-direction (PSI) ! ! JstrT Starting nest tile in the J-direction (RHO-points) ! ! JendT Ending nest tile in the J-direction (RHO-points) ! ! ! ! Istrm3 Starting private I-halo computations, Istr-3 ! ! Istrm2 Starting private I-halo computations, Istr-2 ! ! Istrm1 Starting private I-halo computations, Istr-1 ! ! IstrUm2 Starting private I-halo computations, IstrU-2 ! ! IstrUm1 Starting private I-halo computations, IstrU-1 ! ! Iendp1 Ending private I-halo computations, Iend+1 ! ! Iendp2 Ending private I-halo computations, Iend+2 ! ! Iendp2i Ending private I-halo computations, Iend+2 interior ! ! Iendp3 Ending private I-halo computations, Iend+3 ! ! ! ! Jstrm3 Starting private J-halo computations, Jstr-3 ! ! Jstrm2 Starting private J-halo computations, Jstr-2 ! ! Jstrm1 Starting private J-halo computations, Jstr-1 ! ! JstrVm2 Starting private J-halo computations, JstrV-2 ! ! JstrVm1 Starting private J-halo computations, JstrV-1 ! ! Jendp1 Ending private J-halo computations, Jend+1 ! ! Jendp2 Ending private J-halo computations, Jend+2 ! ! Jendp2i Ending private J-halo computations, Jend+2 interior ! ! Jendp3 Ending private J-halo computations, Jend+3 ! ! ! !======================================================================! ! USE mod_param ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile integer, intent(out) :: Itile, Jtile integer, intent(out) :: Iend, Istr, Jend, Jstr integer, intent(out) :: IstrM, IstrR, IstrU, IendR integer, intent(out) :: JstrM, JstrR, JstrV, JendR integer, intent(out) :: IstrB, IendB, IstrP, IendP, IstrT, IendT integer, intent(out) :: JstrB, JendB, JstrP, JendP, JstrT, JendT integer, intent(out) :: Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1 integer, intent(out) :: Iendp1, Iendp2, Iendp2i, Iendp3 integer, intent(out) :: Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1 integer, intent(out) :: Jendp1, Jendp2, Jendp2i, Jendp3 ! ! Local variable declarations. ! integer :: my_Istr, my_Iend, my_Jstr, my_Jend ! !----------------------------------------------------------------------- ! Set physical non-overlapping grid bounds according to tile rank. !----------------------------------------------------------------------- ! ! Non-tiled grid bounds. This is used in serial or shared-memory ! modes to compute values in the full grid outside of parallel ! regions. ! IF (tile.eq.-1) THEN Itile=-1 Jtile=-1 my_Istr=1 my_Iend=Lm(ng) my_Jstr=1 my_Jend=Mm(ng) ! ! Tiled grids bounds. ! ELSE CALL tile_bounds_2d (ng, tile, Lm(ng), Mm(ng), Itile, Jtile, & & my_Istr, my_Iend, my_Jstr, my_Jend) END IF ! ! Compute C-staggered variables bounds from tile bounds. ! CALL var_bounds (ng, tile, my_Istr, my_Iend, my_Jstr, my_Jend, & & Istr, Iend, Jstr, Jend, & & IstrM, IstrR, IstrU, IendR, & & JstrM, JstrR, JstrV, JendR, & & IstrB, IendB, IstrP, IendP, IstrT, IendT, & & JstrB, JendB, JstrP, JendP, JstrT, JendT, & & Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1, & & Iendp1, Iendp2, Iendp2i, Iendp3, & & Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1, & & Jendp1, Jendp2, Jendp2i, Jendp3) RETURN END SUBROUTINE get_tile SUBROUTINE tile_bounds_1d (ng, tile, Imax, Istr, Iend) ! !======================================================================= ! ! ! This routine computes the starting and ending indices for the 1D ! ! decomposition between all available threads or partitions. ! ! ! ! 1 _____________________ Imax ! ! ! ! On Input: ! ! ! ! ng Nested grid number ! ! tile Thread or partition ! ! Imax Global number of points ! ! ! ! On Output: ! ! ! ! Istr Starting partition index ! ! Iend Ending partition index ! ! ! !======================================================================! ! USE mod_param ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, Imax integer, intent(out) :: Iend, Istr ! ! Local variable declarations. ! integer :: ChunkSize, Margin, Nnodes ! !----------------------------------------------------------------------- ! Compute 1D decomposition starting and ending indices. !----------------------------------------------------------------------- ! Nnodes=NtileI(ng)*NtileJ(ng) ChunkSize=(Imax+Nnodes-1)/Nnodes Margin=(Nnodes*ChunkSize-Imax)/2 IF (Imax.ge.Nnodes) THEN Istr=1+tile*ChunkSize-Margin Iend=Istr+ChunkSize-1 Istr=MAX(Istr,1) Iend=MIN(Iend,Imax) ELSE Istr=1 Iend=Imax END IF RETURN END SUBROUTINE tile_bounds_1d SUBROUTINE tile_bounds_2d (ng, tile, Imax, Jmax, Itile, Jtile, & & Istr, Iend, Jstr, Jend) ! !======================================================================= ! ! ! This routine computes the starting and ending horizontal indices ! ! for each sub-domain partition or tile for a grid bounded between ! ! (1,1) and (Imax,Jmax): ! ! ! ! _________ (Imax,Jmax) ! ! | | ! ! | | ! ! |_________| ! ! (1,1) ! ! ! ! On Input: ! ! ! ! ng Nested grid number ! ! tile Sub-domain partition ! ! Imax Global number of points in the I-direction ! ! Jmax Global number of points in the J-direction ! ! ! ! On Output: ! ! ! ! Itile I-tile coordinate (a value from 0 to NtileI(ng)) ! ! Jtile J-tile coordinate (a value from 0 to NtileJ(ng)) ! ! Istr Starting tile index in the I-direction ! ! Iend Ending tile index in the I-direction ! ! Jstr Starting tile index in the J-direction ! ! Jend Ending tile index in the J-direction ! ! ! !======================================================================! ! USE mod_param ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, Imax, Jmax integer, intent(out) :: Itile, Jtile integer, intent(out) :: Iend, Istr, Jend, Jstr ! ! Local variable declarations. ! integer :: ChunkSizeI, ChunkSizeJ, MarginI, MarginJ ! !----------------------------------------------------------------------- ! Compute tile decomposition for a horizontal grid bounded between ! (1,1) and (Imax,Jmax). !----------------------------------------------------------------------- ! ChunkSizeI=(Imax+NtileI(ng)-1)/NtileI(ng) ChunkSizeJ=(Jmax+NtileJ(ng)-1)/NtileJ(ng) MarginI=(NtileI(ng)*ChunkSizeI-Imax)/2 MarginJ=(NtileJ(ng)*ChunkSizeJ-Jmax)/2 Jtile=tile/NtileI(ng) Itile=tile-Jtile*NtileI(ng) ! ! Tile bounds in the I-direction. ! Istr=1+Itile*ChunkSizeI-MarginI Iend=Istr+ChunkSizeI-1 Istr=MAX(Istr,1) Iend=MIN(Iend,Imax) ! ! Tile bounds in the J-direction. ! Jstr=1+Jtile*ChunkSizeJ-MarginJ Jend=Jstr+ChunkSizeJ-1 Jstr=MAX(Jstr,1) Jend=MIN(Jend,Jmax) RETURN END SUBROUTINE tile_bounds_2d SUBROUTINE var_bounds (ng, tile, & & my_Istr, my_Iend, my_Jstr, my_Jend, & & Istr, Iend, Jstr, Jend, & & IstrM, IstrR, IstrU, IendR, & & JstrM, JstrR, JstrV, JendR, & & IstrB, IendB, IstrP, IendP, IstrT, IendT, & & JstrB, JendB, JstrP, JendP, JstrT, JendT, & & Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1, & & Iendp1, Iendp2, Iendp2i, Iendp3, & & Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1, & & Jendp1, Jendp2, Jendp2i, Jendp3) ! !======================================================================= ! ! ! This routine computes the computational grid starting and ending ! ! horizontal indices for each C-staggered variable in terms of the ! ! physical grid sub-domain partition or tile. ! ! ! ! On Input: ! ! ! ! ng Nested grid number ! ! tile Domain partition ! ! my_Istr Physical grid starting tile index in the I-direction ! ! my_Iend Physical grid ending tile index in the I-direction ! ! my_Jstr Physical grid starting tile index in the J-direction ! ! my_Jend Physical grid ending tile index in the J-direction ! ! ! ! On Output: ! ! ! ! Istr Starting tile index in the I-direction ! ! Iend Ending tile index in the I-direction ! ! Jstr Starting tile index in the J-direction ! ! Jend Ending tile index in the J-direction ! ! ! ! IstrR Starting tile index in the I-direction (RHO-points) ! ! IstrU Starting tile index in the I-direction (U-points) ! ! IendR Ending tile index in the I-direction (RHO_points) ! ! ! ! JstrR Starting tile index in the J-direction (RHO-points) ! ! JstrV Starting tile index in the J-direction (V-points) ! ! JendR Ending tile index in the J-direction (RHO_points) ! ! ! ! IstrB Starting nest tile in the I-direction (RHO-, V-obc) ! ! IendB Ending nest tile in the I-direction (RHO-, V-obc) ! ! IstrM Starting nest tile in the I-direction (PSI-, U-obc) ! ! IstrP Starting nest tile in the I-direction (PSI, U-points) ! ! IendP Ending nest tile in the I-direction (PSI) ! ! IstrT Starting nest tile in the I-direction (RHO-points) ! ! IendT Ending nest tile in the I-direction (RHO_points) ! ! ! ! JstrB Starting nest tile in the J-direction (RHO-, U-obc) ! ! JendB Ending nest tile in the J-direction (RHO-, U-obc) ! ! JstrM Starting nest tile in the J-direction (PSI-, V-obc) ! ! JstrP Starting nest tile in the J-direction (PSI, V-points) ! ! JendP Ending nest tile in the J-direction (PSI) ! ! JstrT Starting nest tile in the J-direction (RHO-points) ! ! JendT Ending nest tile in the J-direction (RHO-points) ! ! ! ! Istrm3 Starting private I-halo computations, Istr-3 ! ! Istrm2 Starting private I-halo computations, Istr-2 ! ! Istrm1 Starting private I-halo computations, Istr-1 ! ! IstrUm2 Starting private I-halo computations, IstrU-2 ! ! IstrUm1 Starting private I-halo computations, IstrU-1 ! ! Iendp1 Ending private I-halo computations, Iend+1 ! ! Iendp2 Ending private I-halo computations, Iend+2 ! ! Iendp2i Ending private I-halo computations, Iend+2 interior ! ! Iendp3 Ending private I-halo computations, Iend+3 ! ! ! ! Jstrm3 Starting private J-halo computations, Jstr-3 ! ! Jstrm2 Starting private J-halo computations, Jstr-2 ! ! Jstrm1 Starting private J-halo computations, Jstr-1 ! ! JstrVm2 Starting private J-halo computations, JstrV-2 ! ! JstrVm1 Starting private J-halo computations, JstrV-1 ! ! Jendp1 Ending private J-halo computations, Jend+1 ! ! Jendp2 Ending private J-halo computations, Jend+2 ! ! Jendp2i Ending private J-halo computations, Jend+2 interior ! ! Jendp3 Ending private J-halo computations, Jend+3 ! ! ! !======================================================================! ! USE mod_param USE mod_scalars ! implicit none ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile integer, intent(in) :: my_Istr, my_Iend, my_Jstr, my_Jend ! integer, intent(out) :: Iend, Istr, Jend, Jstr integer, intent(out) :: IstrM, IstrR, IstrU, IendR integer, intent(out) :: JstrM, JstrR, JstrV, JendR integer, intent(out) :: IstrB, IendB, IstrP, IendP, IstrT, IendT integer, intent(out) :: JstrB, JendB, JstrP, JendP, JstrT, JendT integer, intent(out) :: Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1 integer, intent(out) :: Iendp1, Iendp2, Iendp2i, Iendp3 integer, intent(out) :: Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1 integer, intent(out) :: Jendp1, Jendp2, Jendp2i, Jendp3 ! !======================================================================= ! Compute lower and upper bounds over a particular domain partition or ! tile for RHO-, U-, and V-variables. !======================================================================= ! ! ROMS uses at staggered stencil: ! ! -------v(i,j+1,k)------- ------W(i,j,k)------- ! | | | | ! u(i,j,k) r(i,j,k) u(i+1,j,k) | r(i,j,k) | ! | | | | ! --------v(i,j,k)-------- -----W(i,j,k-1)------ ! ! horizontal stencil vertical stencil ! C-grid ! ! ! M r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r ! : : ! M v p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p v ! : + | | | | | | | | + : ! Mm r u r u r u r u r u r u r u r u r u r u r ! : + | | | | | | | | + : ! Mm v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v ! : + | | | | | | | | + : ! r u r u r u r u r u r u r u r u r u r u r ! : + | | | | | | | | + : ! v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v ! : + | | | | | | | | + : ! r u r u r u r u r u r u r u r u r u r u r ! : + | | | | | | | | + : ! v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v ! : + | | | | | | | | + : ! 2 r u r u r u r u r u r u r u r u r u r u r ! : + | | | | | | | | + : ! 2 v p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p--v--p v ! : + | | | | | | | | + : ! 1 r u r u r u r u r u r u r u r u r u r u r ! : + | | | | | | | | + : ! 1 v p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p++v++p v ! : : ! 0 r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r..u..r ! 1 2 Lm L ! 0 1 2 Lm L ! ! interior Boundary Conditions ! computations W E S N ! ! RH0-type variables: [1:Lm, 1:Mm] [0,:] [L,:] [:,0] [:,M] ! PSI-type variables: [2:Lm, 2:Mm] [1,:] [L,:] [:,1] [:,M] ! U-type variables: [2:Lm, 1:Mm] [1,:] [L,:] [:,0] [:,M] ! V-type variables: [1:Lm, 2:Mm] [0,:] [L,:] [:,1] [:,M] ! ! Compute derived bounds for the loop indices over a subdomain tile. ! The extended bounds (labelled by suffix R) are designed to cover ! also the outer grid points (outlined above with :), if the subdomain ! tile is adjacent to the physical boundary (outlined above with +). ! Notice that IstrR, IendR, JstrR, JendR tile bounds computed here ! DO NOT COVER ghost points (outlined below with *) associated with ! periodic boundaries (if any) or the computational margins of 1 ! subdomains. ! ! Left/Top Tile Right/Top Tile ! ! JendR r..u..r..u..r..u..r..u * * * * u..r..u..r..u..r..u..r ! : Istr Iend Istr Iend : ! v p++v++p++v++p++v++p * * Jend * * p++v++p++v++p++v++p v ! : + | | | | | | + : ! r u r u r u r u * * * * u r u r u r u r ! : + | | | | | | + : ! v p--v--p--v--p--v--p * * * * p--v--p--v--p--v--p v ! : + | | | | | | + : ! r u r u r u r u * * * * u r u r u r u r ! : + | | | | | | + : ! v p--v--p--v--p--v--p * * Jstr * * p--v--p--v--p--v--p v ! ! * * * * * * * * * * * * * * * * * * * * ! ! * * * * * * * * * * * * * * * * * * * * ! ! IstrR IstrU IendR ! IstrM ! ! IstrB=Istr ! * * * * * * * * * * * IstrM=Istr ! Ghost Points IstrP=Istr ! * * * * * * * * * * * IstrR=Istr ! IstrT=Istr ! * * p--v--p--v--p--v--p * * Jend IstrU=Istr ! | | | | IendB=Iend ! Interior * * u r u r u r u * * IendP=Iend ! Tile | | | | IendR=Iend ! * * p--v--p--v--p--v--p * * IendT=Iend ! | | | | JstrB=Jstr ! * * u r u r u r u * * JstrM=Jstr ! | | | | JstrP=Jstr ! * * p--v--p--v--p--v--p * * Jstr JstrR=Jstr ! JstrT=Jstr ! * * * * * * * * * * * JstrV=Jstr ! JendB=Jend ! * * * * * * * * * * * JendP=Jend ! JendR=Jend ! Istr Iend JendT=Jend ! ! ! ! * * * * * * * * * * * * * * * * * * * * ! ! * * * * * * * * * * * * * * * * * * * * ! Istr Iend ! v p--v--p--v--p--v--p * * Jend * * p--v--p--v--p--v--p v ! : + | | | | | | + : ! r u r u r u r u * * * * u r u r u r u r ! : + | | | | | | + : ! JstrV v p--v--p--v--p--v--p * * * * p--v--p--v--p--v--p v JstrM ! : + | | | | | | + : ! r u r u r u r u * * * * u r u r u r u r ! : + | | | | | | + : ! v p++v++p++v++p++v++p * * Jstr * * p++v++p++v++p++v++p v ! : : ! r..u..r..u..r..u..r..u * * * * u..r..u..r..u..r..u..r ! ! IstrR IstrU IendR ! IstrM ! ! Left/Bottom Tile Right/Bottom Tile ! ! ! It also computes loop-bounds for U- and V-type variables which ! belong to the interior of the computational domain. These are ! labelled by suffixes U,V and they step one grid point inward from ! the side of the subdomain adjacent to the physical boundary. ! Conversely, for an internal subdomain which does not include a ! segment of the physical boundary, all bounds with suffixes R,U,V ! are set to the same values of corresponding non-suffixed bounds. ! ! In nested grids there are additional indices to process the ! overlap regions and ranges when calling lateral boundary ! conditions routines (see diagrams below): ! ! IstrT starting overlap I-direction (RHO-points) ! IendT ending overlap I-direction (RHO-points) ! JstrT starting overlap J-direction (RHO-points) ! JendT ending overlap J-direction (RHO-points) ! ! IstrP starting overlap I-direction (PSI-, U-points) ! IendP ending overlap I-direction (PSI-points) ! JstrP starting overlap J-direction (PSI-, V-points) ! JendP ending overlap J-direction (PSI-points) ! ! IstrB starting boundary I-direction (RHO-, V-points), IstrT+1 ! IendB ending boundary I-direction (RHO-, V-points), IendT-1 ! JstrB starting boundary J-direction (RHO-, U-points), JstrT+1 ! JendB ending boundary J-direction (RHO-, U-points), JendT-1 ! ! IstrM starting boundary I-direction (PSI-, U-points), IstrP+1 ! JstrM starting boundary J-direction (PSI-, V-points), JstrP+1 ! ! If not nesting, these indices are set to: ! ! IstrT = IstrR IstrP = Istr IstrB = Istr IstrM = IstrU ! IendT = IendR IendP = Iend IendB = Iend ! JstrT = JstrR JstrP = Jstr JstrB = Jstr JstrM = JstrV ! JendT = JendR JendP = Jend JendB = Jend ! ! The following diagram shows the lower left corner in nested grids: ! ! + ! r u r u r u r u r u r u r u r u r u r ! : : : + | | | | | ! v..p..v..p..v..p..v p--v--p--v--p--v--p--v--p--v--p--v ! : : : + | | | | | ! 2 r u r u r u r u r u r u r u r u r u r ! : : : + | | | | | ! 2 v..p..v..p..v..p..v p--v--p--v--p--v--p--v--p--v--p--v JstrV ! : : : + | | | | | ! 1 r u r u r u r u r u r u r u r u r u r Jstr ! : : : + | | | | | ! 1 v..p..v..p..v..p..v p++v++p++v++p++v++p++v++p++v++p++v++ ! : : : ! 0 r u r u r u r u r u r u r u r u r u r JstrR ! : : : : : : : : : ! 0 v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v ! : : : : : : : : : ! -1 r u r u r u r u r u r u r u r u r u r ! : : : : : : : : : ! -1 v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v JstrM ! : : : : : : : : : ! -2 r u r u r u r u r u r u r u r u r u r JstrB ! : : : : : : : : : ! -2 v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v JstrP ! : : : : : : : : : ! -3 r u r u r u r u r u r u r u r u r u r JstrT ! ! -2 -1 0 1 2 ! -3 -2 -1 0 1 2 ! ! IstrP IstrM IstrU ! IstrB Istr ! IstrT IstrR ! ! ! The following diagram shows the upper right in nested grids: ! ! ! IendR IendT ! Iend IendB ! IendP ! ! Lm L L+1 L+2 ! Lm L L+1 L+2 ! ! M+2 r u r u r u r u r u r u r u r u r JendT ! : : : : : : : : ! M+1 v p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v JendP ! : : : : : : : : ! M+1 r u r u r u r u r u r u r u r u r JendB ! : : : : : : : : ! M+1 v p..v..p..v..p..v..p..v..p..v..p..v..p..v..p..v ! : : : : : : : : ! M r u r u r u r u r u r u r u r u r JendR ! : : ! M ++v++p++v++p++v++p++v++p++v++p++v++p v..p..v..p..v ! | | | | | + : : ! Mm r u r u r u r u r u r u r u r u r Jend ! | | | | | + : : ! Mm v--p--v--p--v--p--v--p--v--p--v--p v..p..v..p..v ! | | | | | + : : ! r u r u r u r u r u r u r u r u r ! | | | | | + : : ! v--p--v--p--v--p--v--p--v--p--v--p v..p..v..p..v ! | | | | | + : : ! r u r u r u r u r u r u r u r u r ! + ! !----------------------------------------------------------------------- ! Starting I-tile indices. !----------------------------------------------------------------------- ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Istr =my_Istr IstrP=my_Istr IstrR=my_Istr IstrT=IstrR IstrU=my_Istr IstrB=my_Istr IstrM=IstrU ELSE Istr =my_Istr IstrP=my_Istr IstrR=my_Istr-1 IstrT=IstrR IstrU=my_Istr+1 IstrB=IstrT+1 IstrM=IstrP+1 END IF ELSE Istr =my_Istr IstrP=my_Istr IstrR=my_Istr IstrT=IstrR IstrU=my_Istr IstrB=my_Istr IstrM=IstrU END IF ! ! Special case, Istrm3: used when MAX(0,Istr-3) is needed. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Istrm3=my_Istr-3 ELSE Istrm3=MAX(0,my_Istr-3) END IF ELSE Istrm3=my_Istr-3 END IF ! ! Special case, Istrm2: used when MAX(0,Istr-2) is needed. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Istrm2=my_Istr-2 ELSE Istrm2=MAX(0,my_Istr-2) END IF ELSE Istrm2=my_Istr-2 END IF ! ! Special case, IstrUm2: used when MAX(1,IstrU-2) is needed. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN IstrUm2=IstrU-2 ELSE IstrUm2=MAX(1,IstrU-2) END IF ELSE IstrUm2=IstrU-2 END IF ! ! Special case, Istrm1: used when MAX(1,Istr-1) is needed. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Istrm1=my_Istr-1 ELSE Istrm1=MAX(1,my_Istr-1) END IF ELSE Istrm1=my_Istr-1 END IF ! ! Special case, IstrUm1: used when MAX(2,IstrU-1) is needed. ! IF (DOMAIN(ng)%Western_Edge(tile)) THEN IF (EWperiodic(ng)) THEN IstrUm1=IstrU-1 ELSE IstrUm1=MAX(2,IstrU-1) END IF ELSE IstrUm1=IstrU-1 END IF ! !----------------------------------------------------------------------- ! Ending I-tile indices. !----------------------------------------------------------------------- ! IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Iend =my_Iend IendR=my_Iend IendP=IendR ! check this one, same as IendR? IendT=IendR IendB=my_Iend ELSE Iend =my_Iend IendR=my_Iend+1 IendP=IendR IendT=IendR IendB=IendT-1 END IF ELSE Iend =my_Iend IendR=my_Iend IendP=IendR IendT=IendR IendB=my_Iend END IF ! ! Special case, Iendp1: used when MIN(Iend+1,Lm(ng)) is needed. ! IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Iendp1=my_Iend+1 ELSE Iendp1=MIN(my_Iend+1,Lm(ng)) END IF ELSE Iendp1=my_Iend+1 END IF ! ! Special case, Iendp2i: used when MIN(Iend+2,Lm(ng)) is needed. ! IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Iendp2i=my_Iend+2 ELSE Iendp2i=MIN(my_Iend+2,Lm(ng)) END IF ELSE Iendp2i=my_Iend+2 END IF ! ! Special case, Iendp2: used when MIN(Iend+2,Lm(ng)+1) is needed. ! IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Iendp2=my_Iend+2 ELSE Iendp2=MIN(my_Iend+2,Lm(ng)+1) END IF ELSE Iendp2=my_Iend+2 END IF ! ! Special case, Iendp3: used when MIN(Iend+3,Lm(ng)+1) is needed. ! IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN IF (EWperiodic(ng)) THEN Iendp3=my_Iend+3 ELSE Iendp3=MIN(my_Iend+3,Lm(ng)+1) END IF ELSE Iendp3=my_Iend+3 END IF ! !----------------------------------------------------------------------- ! Starting J-tile indices. !----------------------------------------------------------------------- ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jstr =my_Jstr JstrP=my_Jstr JstrR=my_Jstr JstrT=JstrR JstrV=my_Jstr JstrB=my_Jstr JstrM=JstrV ELSE Jstr =my_Jstr JstrP=my_Jstr JstrR=my_Jstr-1 JstrT=JstrR JstrV=my_Jstr+1 JstrB=JstrT+1 JstrM=JstrP+1 END IF ELSE Jstr =my_Jstr JstrP=my_Jstr JstrR=my_Jstr JstrT=JstrR JstrV=my_Jstr JstrB=my_Jstr JstrM=JstrV END IF ! ! Special case, Jstrm3: used when MAX(0,Jstr-3) is needed. ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jstrm3=my_Jstr-3 ELSE Jstrm3=MAX(0,my_Jstr-3) END IF ELSE Jstrm3=my_Jstr-3 END IF ! ! Special case, Jstrm2: used when MAX(0,Jstr-2) is needed. ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jstrm2=my_Jstr-2 ELSE Jstrm2=MAX(0,my_Jstr-2) END IF ELSE Jstrm2=my_Jstr-2 END IF ! ! Special case, JstrVm2: used when MAX(1,JstrV-2) is needed. ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN JstrVm2=JstrV-2 ELSE JstrVm2=MAX(1,JstrV-2) END IF ELSE JstrVm2=JstrV-2 END IF ! ! Special case, Jstrm1: used when MAX(1,Jstr-1) is needed. ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jstrm1=my_Jstr-1 ELSE Jstrm1=MAX(1,my_Jstr-1) END IF ELSE Jstrm1=my_Jstr-1 END IF ! ! Special case, JstrVm1: used when MAX(2,JstrV-1) is needed. ! IF (DOMAIN(ng)%Southern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN JstrVm1=JstrV-1 ELSE JstrVm1=MAX(2,JstrV-1) END IF ELSE JstrVm1=JstrV-1 END IF ! !----------------------------------------------------------------------- ! Ending J-tile indices. !----------------------------------------------------------------------- ! IF (DOMAIN(ng)%Northern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jend =my_Jend JendR=my_Jend JendP=JendR ! check this one, same as JendR? JendT=JendR JendB=my_Jend ELSE Jend =my_Jend JendR=my_Jend+1 JendP=JendR JendT=JendR JendB=JendT-1 END IF ELSE Jend =my_Jend JendR=my_Jend JendP=JendR JendP=JendR JendT=JendR JendB=my_Jend END IF ! ! Special case, Jendp1: used when MIN(Jend+1,Mm(ng)) is needed. ! IF (DOMAIN(ng)%Northern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jendp1=my_Jend+1 ELSE Jendp1=MIN(my_Jend+1,Mm(ng)) END IF ELSE Jendp1=my_Jend+1 END IF ! ! Special case, Jendp2i: used when MIN(Jend+2,Mm(ng)) is needed. ! IF (DOMAIN(ng)%Northern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jendp2i=my_Jend+2 ELSE Jendp2i=MIN(my_Jend+2,Mm(ng)) END IF ELSE Jendp2i=my_Jend+2 END IF ! ! Special case, Jendp2: used when MIN(Jend+2,Mm(ng)+1) is needed. ! IF (DOMAIN(ng)%Northern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jendp2=my_Jend+2 ELSE Jendp2=MIN(my_Jend+2,Mm(ng)+1) END IF ELSE Jendp2=my_Jend+2 END IF ! ! Special case, Jendp3: used when MIN(Jend+3,Mm(ng)+1) is needed. ! IF (DOMAIN(ng)%Northern_Edge(tile)) THEN IF (NSperiodic(ng)) THEN Jendp3=my_Jend+3 ELSE Jendp3=MIN(my_Jend+3,Mm(ng)+1) END IF ELSE Jendp3=my_Jend+3 END IF RETURN END SUBROUTINE var_bounds