23 module block_control_mod
27 use fms_string_utils_mod,
only:
string
35 integer,
dimension(:,:),
allocatable :: ix
41 integer,
dimension(:),
allocatable :: ii
42 integer,
dimension(:),
allocatable :: jj
48 integer :: nx_block, ny_block
50 integer :: isc, iec, jsc, jec
52 integer,
dimension(:),
allocatable :: ibs , & !< block extents for mpp-style
56 type(
ix_type),
dimension(:),
allocatable :: ix
58 integer,
dimension(:),
allocatable :: blksz
60 integer,
dimension(:,:),
allocatable :: blkno
61 integer,
dimension(:,:),
allocatable :: ixp
63 type(
pk_type),
dimension(:),
allocatable :: index
79 nx_block, ny_block, message)
80 character(len=*),
intent(in) :: component
82 integer,
intent(in) :: isc, iec, jsc, jec, kpts
83 integer,
intent(in) :: nx_block, ny_block
84 logical,
intent(inout) :: message
103 integer,
dimension(nx_block) :: i1, i2
104 integer,
dimension(ny_block) :: j1, j2
105 character(len=256) :: text
106 integer :: i, j, nblks, ix, ii, jj
107 integer :: non_uniform_blocks
110 non_uniform_blocks = 0
111 if ((mod(iec-isc+1,nx_block) .ne. 0) .or. (mod(jec-jsc+1,ny_block) .ne. 0))
then
112 non_uniform_blocks = 1
114 call mpp_sum(non_uniform_blocks)
115 if (non_uniform_blocks > 0 )
then
117 "have non-uniform blocks for block size ("//
string(nx_block)//
","//
string(ny_block)//
")")
123 if (iec-isc+1 .lt. nx_block) &
124 call mpp_error(fatal,
'block_control: number of '//trim(component)//.gt.
' nxblocks &
125 &number of elements in MPI-domain size')
126 if (jec-jsc+1 .lt. ny_block) &
127 call mpp_error(fatal,
'block_control: number of '//trim(component)//.gt.
' nyblocks &
128 &number of elements in MPI-domain size')
132 nblks = nx_block*ny_block
138 block%nx_block = nx_block
139 block%ny_block = ny_block
142 if (.not.
allocated(block%ibs)) &
143 allocate (block%ibs(nblks), &
153 block%ibs(blocks) = i1(i)
154 block%jbs(blocks) = j1(j)
155 block%ibe(blocks) = i2(i)
156 block%jbe(blocks) = j2(j)
157 allocate(block%ix(blocks)%ix(i1(i):i2(i),j1(j):j2(j)) )
162 block%ix(blocks)%ix(ii,jj) = ix
178 kpts, blksz, message)
179 character(len=*),
intent(in) :: component
181 integer,
intent(in) :: isc, iec, jsc, jec, kpts
182 integer,
intent(inout) :: blksz
183 logical,
intent(inout) :: message
197 integer :: nblks, lblksz, tot_pts, nb, ix, ii, jj
198 character(len=256) :: text
200 tot_pts = (iec - isc + 1) * (jec - jsc + 1)
205 nblks = tot_pts/blksz
206 if (mod(tot_pts,blksz) .gt. 0)
then
212 if (mod(tot_pts,blksz) .ne. 0)
then
213 write( text,
'(a,a,2i4,a,i4,a,i4)' ) trim(component),
'define_blocks_packed: domain (',&
214 (iec-isc+1), (jec-jsc+1),
') is not an even divisor with definition (',&
215 blksz,
') - blocks will not be uniform with a remainder of ',mod(tot_pts,blksz)
227 if (.not.
allocated(block%blksz)) &
228 allocate (block%blksz(nblks), &
229 block%index(nblks), &
230 block%blkno(isc:iec,jsc:jec), &
231 block%ixp(isc:iec,jsc:jec))
236 if (nb .EQ. nblks) lblksz = tot_pts - (nb-1) * blksz
237 block%blksz(nb) = lblksz
238 allocate (block%index(nb)%ii(lblksz), &
239 block%index(nb)%jj(lblksz))
248 if (ix .GT. blksz)
then
252 block%ixp(ii,jj) = ix
253 block%blkno(ii,jj) = nb
254 block%index(nb)%ii(ix) = ii
255 block%index(nb)%jj(ix) = jj
261 end module block_control_mod
subroutine, public define_blocks_packed(component, Block, isc, iec, jsc, jec, kpts, blksz, message)
Creates and populates a data type which is used for defining the sub-blocks of the MPI-domain to enha...
subroutine, public define_blocks(component, Block, isc, iec, jsc, jec, kpts, nx_block, ny_block, message)
Sets up "blocks" used for OpenMP threading of column-based calculations using rad_n[x/y]xblock from c...
Block data and extents for OpenMP threading of column-based calculations.
Type to dereference packed index from global index.
Type to dereference packed index from global indices.
character(:) function, allocatable, public string(v, fmt)
Converts a number or a Boolean value to a string.
subroutine mpp_compute_extent(isg, ieg, ndivs, ibegin, iend, extent)
Computes extents for a grid decomposition with the given indices and divisions.
integer function mpp_npes()
Returns processor count for current pelist.