FMS  2025.04
Flexible Modeling System
mpp_domains_util.inc
1 ! -*-f90-*-
2 
3 !***********************************************************************
4 !* Apache License 2.0
5 !*
6 !* This file is part of the GFDL Flexible Modeling System (FMS).
7 !*
8 !* Licensed under the Apache License, Version 2.0 (the "License");
9 !* you may not use this file except in compliance with the License.
10 !* You may obtain a copy of the License at
11 !*
12 !* http://www.apache.org/licenses/LICENSE-2.0
13 !*
14 !* FMS is distributed in the hope that it will be useful, but WITHOUT
15 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
16 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
17 !* PARTICULAR PURPOSE. See the License for the specific language
18 !* governing permissions and limitations under the License.
19 !***********************************************************************
20 
21 !> @file
22 !> @brief Utility routines for getting and setting values in @ref mpp_domains_mod
23 
24 !> @ingroup mpp_domains_mod
25 !> @{
26 
27  !> @brief Set user stack size.
28  !!
29  !> This sets the size of an array that is used for internal storage by
30  !! <TT>mpp_domains</TT>. This array is used, for instance, to buffer the
31  !! data sent and received in halo updates.<br>
32  !! This call has implied global synchronization. It should be
33  !! placed somewhere where all PEs can call it.
35  !set the mpp_domains_stack variable to be at least n LONG words long
36  integer, intent(in) :: n
37  character(len=8) :: text
38 
39  if( n.LE.mpp_domains_stack_size )return
40  if( allocated(mpp_domains_stack) )deallocate(mpp_domains_stack)
41  allocate( mpp_domains_stack(n) )
42  if( allocated(mpp_domains_stack_nonblock) )deallocate(mpp_domains_stack_nonblock)
43  allocate( mpp_domains_stack_nonblock(n) )
44 
45  mpp_domains_stack_size = n
46  write( text,'(i8)' )n
47  if( mpp_pe().EQ.mpp_root_pe() )call mpp_error( note, 'MPP_DOMAINS_SET_STACK_SIZE: stack size set to '//text//'.' )
48 
49  return
50  end subroutine mpp_domains_set_stack_size
51 
52 
53 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
54  ! !
55  ! MPP_DOMAINS: overloaded operators (==, /=) !
56  ! !
57 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
58 
59  function mpp_domain1d_eq( a, b )
60  logical :: mpp_domain1d_eq
61  type(domain1d), intent(in) :: a, b
62 
63  mpp_domain1d_eq = ( a%compute%begin.EQ.b%compute%begin .AND. &
64  a%compute%end .EQ.b%compute%end .AND. &
65  a%domain_data%begin .EQ.b%domain_data%begin .AND. &
66  a%domain_data%end .EQ.b%domain_data%end .AND. &
67  a%global%begin .EQ.b%global%begin .AND. &
68  a%global%end .EQ.b%global%end )
69  !compare pelists
70  ! if( mpp_domain1D_eq )mpp_domain1D_eq = ASSOCIATED(a%list) .AND. ASSOCIATED(b%list)
71  ! if( mpp_domain1D_eq )mpp_domain1D_eq = size(a%list(:)).EQ.size(b%list(:))
72  ! if( mpp_domain1D_eq )mpp_domain1D_eq = ALL(a%list%pe.EQ.b%list%pe)
73 
74  return
75  end function mpp_domain1d_eq
76 
77  function mpp_domain1d_ne( a, b )
78  logical :: mpp_domain1d_ne
79  type(domain1d), intent(in) :: a, b
80 
81  mpp_domain1d_ne = .NOT. ( a.EQ.b )
82  return
83  end function mpp_domain1d_ne
84 
85  function mpp_domain2d_eq( a, b )
86  logical :: mpp_domain2d_eq
87  type(domain2d), intent(in) :: a, b
88  integer :: nt, n
89 
90  mpp_domain2d_eq = size(a%x(:)) .EQ. size(b%x(:))
91  nt = size(a%x(:))
92  do n = 1, nt
93  if(mpp_domain2d_eq) mpp_domain2d_eq = a%x(n).EQ.b%x(n) .AND. a%y(n).EQ.b%y(n)
94  end do
95 
96  if( mpp_domain2d_eq .AND. ((a%pe.EQ.null_pe).OR.(b%pe.EQ.null_pe)) )return !NULL_DOMAIN2D
97  !compare pelists
98  if( mpp_domain2d_eq )mpp_domain2d_eq = ASSOCIATED(a%list) .AND. ASSOCIATED(b%list)
99  if( mpp_domain2d_eq )mpp_domain2d_eq = size(a%list(:)).EQ.size(b%list(:))
100  if( mpp_domain2d_eq )mpp_domain2d_eq = all(a%list%pe.EQ.b%list%pe)
101  if( mpp_domain2d_eq )mpp_domain2d_eq = all(a%io_layout .EQ. b%io_layout)
102  if( mpp_domain2d_eq )mpp_domain2d_eq = a%symmetry .eqv. b%symmetry
103 
104  return
105  end function mpp_domain2d_eq
106 
107  !#####################################################################
108 
109  function mpp_domain2d_ne( a, b )
110  logical :: mpp_domain2d_ne
111  type(domain2d), intent(in) :: a, b
112 
113  mpp_domain2d_ne = .NOT. ( a.EQ.b )
114  return
115  end function mpp_domain2d_ne
116 
117 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
118  ! !
119  ! MPP_GET and SET routiness: retrieve various components of domains !
120  ! !
121 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
122 
123  subroutine mpp_get_compute_domain1d( domain, begin, end, size, max_size, is_global )
124  type(domain1d), intent(in) :: domain
125  integer, intent(out), optional :: begin, end, size, max_size
126  logical, intent(out), optional :: is_global
127 
128  if( PRESENT(begin) )begin = domain%compute%begin
129  if( PRESENT(end) )end = domain%compute%end
130  if( PRESENT(size) )size = domain%compute%size
131  if( PRESENT(max_size) )max_size = domain%compute%max_size
132  if( PRESENT(is_global) )is_global = domain%compute%is_global
133  return
134  end subroutine mpp_get_compute_domain1d
135 
136  !#####################################################################
137  subroutine mpp_get_data_domain1d( domain, begin, end, size, max_size, is_global )
138  type(domain1d), intent(in) :: domain
139  integer, intent(out), optional :: begin, end, size, max_size
140  logical, intent(out), optional :: is_global
141 
142  if( PRESENT(begin) )begin = domain%domain_data%begin
143  if( PRESENT(end) )end = domain%domain_data%end
144  if( PRESENT(size) )size = domain%domain_data%size
145  if( PRESENT(max_size) )max_size = domain%domain_data%max_size
146  if( PRESENT(is_global) )is_global = domain%domain_data%is_global
147  return
148  end subroutine mpp_get_data_domain1d
149 
150  !#####################################################################
151  subroutine mpp_get_global_domain1d( domain, begin, end, size, max_size )
152  type(domain1d), intent(in) :: domain
153  integer, intent(out), optional :: begin, end, size, max_size
154 
155  if( PRESENT(begin) )begin = domain%global%begin
156  if( PRESENT(end) )end = domain%global%end
157  if( PRESENT(size) )size = domain%global%size
158  if( PRESENT(max_size) )max_size = domain%global%max_size
159  return
160  end subroutine mpp_get_global_domain1d
161 
162  !#####################################################################
163  subroutine mpp_get_memory_domain1d( domain, begin, end, size, max_size, is_global )
164  type(domain1d), intent(in) :: domain
165  integer, intent(out), optional :: begin, end, size, max_size
166  logical, intent(out), optional :: is_global
167 
168  if( PRESENT(begin) )begin = domain%memory%begin
169  if( PRESENT(end) )end = domain%memory%end
170  if( PRESENT(size) )size = domain%memory%size
171  if( PRESENT(max_size) )max_size = domain%memory%max_size
172  if( PRESENT(is_global) )is_global = domain%memory%is_global
173  return
174  end subroutine mpp_get_memory_domain1d
175 
176  !#####################################################################
177  subroutine mpp_get_compute_domain2d( domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, &
178  x_is_global, y_is_global, tile_count, position )
179  type(domain2d), intent(in) :: domain
180  integer, intent(out), optional :: xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size
181  logical, intent(out), optional :: x_is_global, y_is_global
182  integer, intent(in), optional :: tile_count, position
183  integer :: tile, ishift, jshift
184 
185  tile = 1
186  if(present(tile_count)) tile = tile_count
187 
188  call mpp_get_compute_domain( domain%x(tile), xbegin, xend, xsize, xmax_size, x_is_global )
189  call mpp_get_compute_domain( domain%y(tile), ybegin, yend, ysize, ymax_size, y_is_global )
190  call mpp_get_domain_shift( domain, ishift, jshift, position )
191  if( PRESENT(xend) ) xend = xend + ishift
192  if( PRESENT(yend) ) yend = yend + jshift
193  if( PRESENT(xsize)) xsize = xsize + ishift
194  if( PRESENT(ysize)) ysize = ysize + jshift
195  if(PRESENT(xmax_size))xmax_size = xmax_size + ishift
196  if(PRESENT(ymax_size))ymax_size = ymax_size + jshift
197 
198  return
199  end subroutine mpp_get_compute_domain2d
200 
201  !#####################################################################
202  subroutine mpp_get_data_domain2d( domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, &
203  x_is_global, y_is_global, tile_count, position )
204  type(domain2d), intent(in) :: domain
205  integer, intent(out), optional :: xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size
206  logical, intent(out), optional :: x_is_global, y_is_global
207  integer, intent(in), optional :: tile_count, position
208  integer :: tile, ishift, jshift
209 
210  tile = 1
211  if(present(tile_count)) tile = tile_count
212 
213  call mpp_get_data_domain( domain%x(tile), xbegin, xend, xsize, xmax_size, x_is_global )
214  call mpp_get_data_domain( domain%y(tile), ybegin, yend, ysize, ymax_size, y_is_global )
215  call mpp_get_domain_shift( domain, ishift, jshift, position )
216  if( PRESENT(xend) ) xend = xend + ishift
217  if( PRESENT(yend) ) yend = yend + jshift
218  if( PRESENT(xsize)) xsize = xsize + ishift
219  if( PRESENT(ysize)) ysize = ysize + jshift
220  if(PRESENT(xmax_size))xmax_size = xmax_size + ishift
221  if(PRESENT(ymax_size))ymax_size = ymax_size + jshift
222 
223  return
224  end subroutine mpp_get_data_domain2d
225 
226  !#####################################################################
227  subroutine mpp_get_global_domain2d( domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, &
228  tile_count, position )
229  type(domain2d), intent(in) :: domain
230  integer, intent(out), optional :: xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size
231  integer, intent(in), optional :: tile_count, position
232  integer :: tile, ishift, jshift
233 
234  tile = 1
235  if(present(tile_count)) tile = tile_count
236 
237  call mpp_get_global_domain( domain%x(tile), xbegin, xend, xsize, xmax_size )
238  call mpp_get_global_domain( domain%y(tile), ybegin, yend, ysize, ymax_size )
239  call mpp_get_domain_shift( domain, ishift, jshift, position )
240  if( PRESENT(xend) ) xend = xend + ishift
241  if( PRESENT(yend) ) yend = yend + jshift
242  if( PRESENT(xsize)) xsize = xsize + ishift
243  if( PRESENT(ysize)) ysize = ysize + jshift
244  if(PRESENT(xmax_size))xmax_size = xmax_size + ishift
245  if(PRESENT(ymax_size))ymax_size = ymax_size + jshift
246 
247  return
248  end subroutine mpp_get_global_domain2d
249 
250  !#####################################################################
251  subroutine mpp_get_memory_domain2d( domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, &
252  x_is_global, y_is_global, position)
253  type(domain2d), intent(in) :: domain
254  integer, intent(out), optional :: xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size
255  logical, intent(out), optional :: x_is_global, y_is_global
256  integer, intent(in), optional :: position
257  integer :: tile, ishift, jshift
258 
259  tile = 1
260 
261  call mpp_get_memory_domain( domain%x(tile), xbegin, xend, xsize, xmax_size, x_is_global )
262  call mpp_get_memory_domain( domain%y(tile), ybegin, yend, ysize, ymax_size, y_is_global )
263  call mpp_get_domain_shift( domain, ishift, jshift, position )
264  if( PRESENT(xend) ) xend = xend + ishift
265  if( PRESENT(yend) ) yend = yend + jshift
266  if( PRESENT(xsize)) xsize = xsize + ishift
267  if( PRESENT(ysize)) ysize = ysize + jshift
268  if(PRESENT(xmax_size))xmax_size = xmax_size + ishift
269  if(PRESENT(ymax_size))ymax_size = ymax_size + jshift
270 
271  return
272  end subroutine mpp_get_memory_domain2d
273 
274  !> @brief Modifies the indices in the domain_axis_spec type to those of the supergrid
275  subroutine mpp_set_super_grid_indices(grid)
276  type(domain_axis_spec), intent(inout) :: grid !< domain_axis_spec type
277 
278  grid%begin = 2*grid%begin-1
279  grid%end = 2*grid%end+1
280  grid%size = grid%end-grid%begin+1
281 
282  end subroutine mpp_set_super_grid_indices
283 
284  !> @brief Modifies the indices of the input domain to create the supergrid domain
285  !!
286  !> This is an example of how to use mpp_create_super_grid_domain
287  !! @code{.F90}
288  !! call mpp_copy_domain(domain_in, domain_out)
289  !! call super_grid_domain(domain_out)
290  !! @endcode
291  !! domain_in is the original domain, domain_out is the domain with the supergrid indices.
292  subroutine mpp_create_super_grid_domain(domain)
293  type(domain2d), intent(inout) :: domain !< Input domain
294 
295  integer :: xbegin !< Begining x indices
296  integer :: ybegin !< Begining y indices
297  integer :: xend !< Ending x indices
298  integer :: yend !< Ending y indices
299  integer :: xsize !< Size of the x domain
300  integer :: ysize !< Size of the y domain
301  integer :: i !< For loops
302 
303  call mpp_get_global_domain(domain, xbegin=xbegin, xend=xend, ybegin=ybegin, yend=yend, xsize=xsize, ysize=ysize)
304  call mpp_set_global_domain (domain, 2*xbegin-1, 2*xend+1, 2*ybegin-1, 2*yend+1, 2*(xend-xbegin)+3, &
305  & 2*(yend-ybegin)+3)
306 
307  call mpp_get_compute_domain(domain, xbegin=xbegin, xend=xend, ybegin=ybegin, yend=yend, xsize=xsize, ysize=ysize)
308  call mpp_set_compute_domain (domain, 2*xbegin-1, 2*xend+1, 2*ybegin-1, 2*yend+1, 2*(xend-xbegin)+3, &
309  & 2*(yend-ybegin)+3)
310 
311  call mpp_get_data_domain(domain, xbegin=xbegin, xend=xend, ybegin=ybegin, yend=yend, xsize=xsize, ysize=ysize)
312  call mpp_set_data_domain (domain, 2*xbegin-1, 2*xend+1, 2*ybegin-1, 2*yend+1, 2*(xend-xbegin)+3, 2*(yend-ybegin)+3)
313 
314  do i=1, size(domain%list(:))
315  call mpp_set_super_grid_indices(domain%list(i-1)%x(1)%global)
316  call mpp_set_super_grid_indices(domain%list(i-1)%y(1)%global)
317 
318  call mpp_set_super_grid_indices(domain%list(i-1)%x(1)%compute)
319  call mpp_set_super_grid_indices(domain%list(i-1)%y(1)%compute)
320 
321  !< There is no data domain in domain%list
322  !call mpp_set_super_grid_indices(domain%list(i-1)%x(1)%domain_data)
323  !call mpp_set_super_grid_indices(domain%list(i-1)%y(1)%domain_data)
324  enddo
325 
326  do i=1, size(domain%x(1)%list)
327  call mpp_set_super_grid_indices(domain%x(1)%list(i-1)%compute)
328  enddo
329 
330  do i=1, size(domain%y(1)%list)
331  call mpp_set_super_grid_indices(domain%y(1)%list(i-1)%compute)
332  enddo
333 
334  end subroutine mpp_create_super_grid_domain
335 
336  !#####################################################################
337  subroutine mpp_set_compute_domain1d( domain, begin, end, size, is_global )
338  type(domain1d), intent(inout) :: domain
339  integer, intent(in), optional :: begin, end, size
340  logical, intent(in), optional :: is_global
341 
342  if(present(begin)) domain%compute%begin = begin
343  if(present(end)) domain%compute%end = end
344  if(present(size)) domain%compute%size = size
345  if(present(is_global)) domain%compute%is_global = is_global
346 
347  end subroutine mpp_set_compute_domain1d
348 
349  !#####################################################################
350  subroutine mpp_set_compute_domain2d( domain, xbegin, xend, ybegin, yend, xsize, ysize, &
351  x_is_global, y_is_global, tile_count )
352  type(domain2d), intent(inout) :: domain
353  integer, intent(in), optional :: xbegin, xend, ybegin, yend, xsize, ysize
354  logical, intent(in), optional :: x_is_global, y_is_global
355  integer, intent(in), optional :: tile_count
356  integer :: tile
357 
358  tile = 1
359  if(present(tile_count)) tile = tile_count
360 
361  call mpp_set_compute_domain(domain%x(tile), xbegin, xend, xsize, x_is_global)
362  call mpp_set_compute_domain(domain%y(tile), ybegin, yend, ysize, y_is_global)
363 
364  end subroutine mpp_set_compute_domain2d
365 
366  !#####################################################################
367  subroutine mpp_set_data_domain1d( domain, begin, end, size, is_global )
368  type(domain1d), intent(inout) :: domain
369  integer, intent(in), optional :: begin, end, size
370  logical, intent(in), optional :: is_global
371 
372  if(present(begin)) domain%domain_data%begin = begin
373  if(present(end)) domain%domain_data%end = end
374  if(present(size)) domain%domain_data%size = size
375  if(present(is_global)) domain%domain_data%is_global = is_global
376 
377  end subroutine mpp_set_data_domain1d
378 
379  !#####################################################################
380  subroutine mpp_set_data_domain2d( domain, xbegin, xend, ybegin, yend, xsize, ysize, &
381  x_is_global, y_is_global, tile_count )
382  type(domain2d), intent(inout) :: domain
383  integer, intent(in), optional :: xbegin, xend, ybegin, yend, xsize, ysize
384  logical, intent(in), optional :: x_is_global, y_is_global
385  integer, intent(in), optional :: tile_count
386  integer :: tile
387 
388  tile = 1
389  if(present(tile_count)) tile = tile_count
390 
391  call mpp_set_data_domain(domain%x(tile), xbegin, xend, xsize, x_is_global)
392  call mpp_set_data_domain(domain%y(tile), ybegin, yend, ysize, y_is_global)
393 
394  end subroutine mpp_set_data_domain2d
395 
396  !#####################################################################
397  subroutine mpp_set_global_domain1d( domain, begin, end, size)
398  type(domain1d), intent(inout) :: domain
399  integer, intent(in), optional :: begin, end, size
400 
401  if(present(begin)) domain%global%begin = begin
402  if(present(end)) domain%global%end = end
403  if(present(size)) domain%global%size = size
404 
405  end subroutine mpp_set_global_domain1d
406 
407  !#####################################################################
408  subroutine mpp_set_global_domain2d( domain, xbegin, xend, ybegin, yend, xsize, ysize, tile_count )
409  type(domain2d), intent(inout) :: domain
410  integer, intent(in), optional :: xbegin, xend, ybegin, yend, xsize, ysize
411  integer, intent(in), optional :: tile_count
412  integer :: tile
413 
414  tile = 1
415  if(present(tile_count)) tile = tile_count
416  call mpp_set_global_domain(domain%x(tile), xbegin, xend, xsize)
417  call mpp_set_global_domain(domain%y(tile), ybegin, yend, ysize)
418 
419  end subroutine mpp_set_global_domain2d
420 
421  !#####################################################################
422 
423  !> @brief Retrieve 1D components of 2D decomposition.
424  !!
425  !! It is sometime necessary to have direct recourse to the domain1D types
426  !! that compose a domain2D object. This call retrieves them.
427  !! @code{.F90}
428  !! call mpp_get_domain_components( domain, x, y )
429  !! @endcode
430  subroutine mpp_get_domain_components( domain, x, y, tile_count )
431  type(domain2d), intent(in) :: domain
432  type(domain1d), intent(inout), optional :: x, y
433  integer, intent(in), optional :: tile_count
434  integer :: tile
435 
436  tile = 1
437  if(present(tile_count)) tile = tile_count
438  if( PRESENT(x) )x = domain%x(tile)
439  if( PRESENT(y) )y = domain%y(tile)
440  return
441  end subroutine mpp_get_domain_components
442 
443  !#####################################################################
444  subroutine mpp_get_compute_domains1d( domain, begin, end, size )
445  type(domain1d), intent(in) :: domain
446  integer, intent(out), optional, dimension(:) :: begin, end, size
447 
448  if( .NOT.module_is_initialized ) &
449  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: must first call mpp_domains_init.' )
450  !we use shape instead of size for error checks because size is used as an argument
451  if( PRESENT(begin) )then
452  if( any(shape(begin).NE.shape(domain%list)) ) &
453  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: begin array size does not match domain.' )
454  begin(:) = domain%list(:)%compute%begin
455  end if
456  if( PRESENT(end) )then
457  if( any(shape(end).NE.shape(domain%list)) ) &
458  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: end array size does not match domain.' )
459  end(:) = domain%list(:)%compute%end
460  end if
461  if( PRESENT(size) )then
462  if( any(shape(size).NE.shape(domain%list)) ) &
463  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: size array size does not match domain.' )
464  size(:) = domain%list(:)%compute%size
465  end if
466  return
467 end subroutine mpp_get_compute_domains1d
468 
469 !#####################################################################
470 subroutine mpp_get_compute_domains2d( domain, xbegin, xend, xsize, ybegin, yend, ysize, position )
471  type(domain2d), intent(in) :: domain
472  integer, intent(out), optional, dimension(:) :: xbegin, xend, xsize, ybegin, yend, ysize
473  integer, intent(in ), optional :: position
474 
475  integer :: i, ishift, jshift
476 
477  call mpp_get_domain_shift( domain, ishift, jshift, position )
478 
479 
480  if( .NOT.module_is_initialized ) &
481  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: must first call mpp_domains_init.' )
482 
483  if( PRESENT(xbegin) )then
484  if( size(xbegin(:)).NE.size(domain%list(:)) ) &
485  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: xbegin array size does not match domain.' )
486  do i = 1, size(xbegin(:))
487  xbegin(i) = domain%list(i-1)%x(1)%compute%begin
488  end do
489  end if
490  if( PRESENT(xend) )then
491  if( size(xend(:)).NE.size(domain%list(:)) ) &
492  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: xend array size does not match domain.' )
493  do i = 1, size(xend(:))
494  xend(i) = domain%list(i-1)%x(1)%compute%end + ishift
495  end do
496  end if
497  if( PRESENT(xsize) )then
498  if( size(xsize(:)).NE.size(domain%list(:)) ) &
499  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: xsize array size does not match domain.' )
500  do i = 1, size(xsize(:))
501  xsize(i) = domain%list(i-1)%x(1)%compute%size + ishift
502  end do
503  end if
504  if( PRESENT(ybegin) )then
505  if( size(ybegin(:)).NE.size(domain%list(:)) ) &
506  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: ybegin array size does not match domain.' )
507  do i = 1, size(ybegin(:))
508  ybegin(i) = domain%list(i-1)%y(1)%compute%begin
509  end do
510  end if
511  if( PRESENT(yend) )then
512  if( size(yend(:)).NE.size(domain%list(:)) ) &
513  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: yend array size does not match domain.' )
514  do i = 1, size(yend(:))
515  yend(i) = domain%list(i-1)%y(1)%compute%end + jshift
516  end do
517  end if
518  if( PRESENT(ysize) )then
519  if( size(ysize(:)).NE.size(domain%list(:)) ) &
520  call mpp_error( fatal, 'MPP_GET_COMPUTE_DOMAINS: ysize array size does not match domain.' )
521  do i = 1, size(ysize(:))
522  ysize(i) = domain%list(i-1)%y(1)%compute%size + jshift
523  end do
524  end if
525  return
526 end subroutine mpp_get_compute_domains2d
527 
528  !#####################################################################
529  subroutine mpp_get_global_domains1d( domain, begin, end, size )
530  type(domain1d), intent(in) :: domain
531  integer, intent(out), optional, dimension(:) :: begin, end, size
532 
533  if( .NOT.module_is_initialized ) &
534  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: must first call mpp_domains_init.' )
535  !we use shape instead of size for error checks because size is used as an argument
536  if( PRESENT(begin) )then
537  if( any(shape(begin).NE.shape(domain%list)) ) &
538  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: begin array size does not match domain.' )
539  begin(:) = domain%list(:)%global%begin
540  end if
541  if( PRESENT(end) )then
542  if( any(shape(end).NE.shape(domain%list)) ) &
543  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: end array size does not match domain.' )
544  end(:) = domain%list(:)%global%end
545  end if
546  if( PRESENT(size) )then
547  if( any(shape(size).NE.shape(domain%list)) ) &
548  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: size array size does not match domain.' )
549  size(:) = domain%list(:)%global%size
550  end if
551  return
552  end subroutine mpp_get_global_domains1d
553 
554 
555 !#####################################################################
556 subroutine mpp_get_global_domains2d( domain, xbegin, xend, xsize, ybegin, yend, ysize, position )
557  type(domain2d), intent(in) :: domain
558  integer, intent(out), optional, dimension(:) :: xbegin, xend, xsize, ybegin, yend, ysize
559  integer, intent(in ), optional :: position
560 
561  integer :: i, ishift, jshift
562 
563  call mpp_get_domain_shift( domain, ishift, jshift, position )
564 
565 
566  if( .NOT.module_is_initialized ) &
567  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: must first call mpp_domains_init.' )
568 
569  if( PRESENT(xbegin) )then
570  if( size(xbegin(:)).NE.size(domain%list(:)) ) &
571  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: xbegin array size does not match domain.' )
572  do i = 1, size(xbegin(:))
573  xbegin(i) = domain%list(i-1)%x(1)%global%begin
574  end do
575  end if
576  if( PRESENT(xend) )then
577  if( size(xend(:)).NE.size(domain%list(:)) ) &
578  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: xend array size does not match domain.' )
579  do i = 1, size(xend(:))
580  xend(i) = domain%list(i-1)%x(1)%global%end + ishift
581  end do
582  end if
583  if( PRESENT(xsize) )then
584  if( size(xsize(:)).NE.size(domain%list(:)) ) &
585  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: xsize array size does not match domain.' )
586  do i = 1, size(xsize(:))
587  xsize(i) = domain%list(i-1)%x(1)%global%size + ishift
588  end do
589  end if
590  if( PRESENT(ybegin) )then
591  if( size(ybegin(:)).NE.size(domain%list(:)) ) &
592  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: ybegin array size does not match domain.' )
593  do i = 1, size(ybegin(:))
594  ybegin(i) = domain%list(i-1)%y(1)%global%begin
595  end do
596  end if
597  if( PRESENT(yend) )then
598  if( size(yend(:)).NE.size(domain%list(:)) ) &
599  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: yend array size does not match domain.' )
600  do i = 1, size(yend(:))
601  yend(i) = domain%list(i-1)%y(1)%global%end + jshift
602  end do
603  end if
604  if( PRESENT(ysize) )then
605  if( size(ysize(:)).NE.size(domain%list(:)) ) &
606  call mpp_error( fatal, 'MPP_GET_GLOBAL_DOMAINS: ysize array size does not match domain.' )
607  do i = 1, size(ysize(:))
608  ysize(i) = domain%list(i-1)%y(1)%global%size + jshift
609  end do
610  end if
611  return
612 end subroutine mpp_get_global_domains2d
613 
614 
615 !#####################################################################
616 subroutine mpp_get_domain_extents1d(domain, xextent, yextent)
617  type(domain2d), intent(in) :: domain
618  integer, dimension(0:), intent(inout) :: xextent, yextent
619  integer :: n
620 
621  if(domain%ntiles .NE. 1) call mpp_error(fatal,"mpp_domains_util.inc(mpp_get_domain_extents1D): "// &
622  "ntiles is more than 1, please use mpp_get_domain_extents2D")
623  if(size(xextent) .NE. size(domain%x(1)%list(:))) call mpp_error(fatal, &
624  & "mpp_domains_util.inc(mpp_get_domain_extents1D): "// &
625  & "size(xextent) does not equal to size(domain%x(1)%list(:)))")
626  if(size(yextent) .NE. size(domain%y(1)%list(:))) call mpp_error(fatal, &
627  & "mpp_domains_util.inc(mpp_get_domain_extents1D): "// &
628  & "size(yextent) does not equal to size(domain%y(1)%list(:)))")
629  do n = 0, size(domain%x(1)%list(:))-1
630  xextent(n) = domain%x(1)%list(n)%compute%size
631  enddo
632  do n = 0, size(domain%y(1)%list(:))-1
633  yextent(n) = domain%y(1)%list(n)%compute%size
634  enddo
635 
636 end subroutine mpp_get_domain_extents1d
637 
638 !#####################################################################
639 !> This will return xextent and yextent for each tile
640 subroutine mpp_get_domain_extents2d(domain, xextent, yextent)
641  type(domain2d), intent(in) :: domain
642  integer, dimension(:,:), intent(inout) :: xextent, yextent
643  integer :: ntile, nlist, n, m, ndivx, ndivy, tile, pos
644 
645  ntile = domain%ntiles
646  nlist = size(domain%list(:))
647  if(size(xextent,2) .ne. ntile .or. size(yextent,2) .ne. ntile) call mpp_error(fatal, &
648  "mpp_domains_utile.inc: the second dimension size of xextent/yextent is not correct")
649  ndivx = size(xextent,1); ndivy = size(yextent,1)
650  do n = 0, nlist-1
651  if(any(domain%list(n)%x(:)%pos>ndivx-1) ) call mpp_error(fatal, &
652  "mpp_domains_utile.inc: first dimension size of xextent is less than the x-layout in some tile")
653  if(any(domain%list(n)%y(:)%pos>ndivy-1) ) call mpp_error(fatal, &
654  "mpp_domains_utile.inc: first dimension size of yextent is less than the y-layout in some tile")
655  end do
656 
657  xextent = 0; yextent=0
658 
659  do n = 0, nlist-1
660  do m = 1, size(domain%list(n)%tile_id(:))
661  tile = domain%list(n)%tile_id(m)
662  pos = domain%list(n)%x(m)%pos+1
663  if(xextent(pos, tile) == 0) xextent(pos,tile) = domain%list(n)%x(m)%compute%size
664  pos = domain%list(n)%y(m)%pos+1
665  if(yextent(pos, tile) == 0) yextent(pos,tile) = domain%list(n)%y(m)%compute%size
666  end do
667  end do
668 
669 
670 end subroutine mpp_get_domain_extents2d
671 
672 !#####################################################################
673 function mpp_get_domain_pe(domain)
674  type(domain2d), intent(in) :: domain
675  integer :: mpp_get_domain_pe
676 
677  mpp_get_domain_pe = domain%pe
678 
679 
680 end function mpp_get_domain_pe
681 
682 
684  type(domain2d), intent(in) :: domain
685  integer :: mpp_get_domain_tile_root_pe
686 
687  mpp_get_domain_tile_root_pe = domain%tile_root_pe
688 
689 end function mpp_get_domain_tile_root_pe
690 
691 
693  type(domain2d), intent(in) :: domain !> domain you are querying for information
694  integer :: mpp_get_domain_tile_commid !> declaration of the return tile communicator
695 
696  mpp_get_domain_tile_commid = domain%tile_comm_id
697 
698 end function mpp_get_domain_tile_commid
699 
700 
701 function mpp_get_domain_commid(domain)
702  type(domain2d), intent(in) :: domain !> domain you are querying for information
703  integer :: mpp_get_domain_commid !> declaration of the return domain communicator
704 
705  mpp_get_domain_commid = domain%comm_id
706 
707 end function mpp_get_domain_commid
708 
709 
710 function mpp_get_io_domain(domain)
711  type(domain2d), intent(in) :: domain
712  type(domain2d), pointer :: mpp_get_io_domain
713 
714  if(ASSOCIATED(domain%io_domain)) then
715  mpp_get_io_domain => domain%io_domain
716  else
717  mpp_get_io_domain => null()
718  endif
719 
720 end function mpp_get_io_domain
721 
722 !#####################################################################
723 ! <SUBROUTINE NAME="mpp_get_pelist1D" INTERFACE="mpp_get_pelist">
724 ! <IN NAME="domain" TYPE="type(domain1D)"></IN>
725 ! <OUT NAME="pelist" TYPE="integer" DIM="(:)"></OUT>
726 ! <OUT NAME="pos" TYPE="integer"></OUT>
727 ! </SUBROUTINE>
728 subroutine mpp_get_pelist1d( domain, pelist, pos )
729  type(domain1d), intent(in) :: domain
730  integer, intent(out) :: pelist(:)
731  integer, intent(out), optional :: pos
732  integer :: ndivs
733 
734  if( .NOT.module_is_initialized ) &
735  call mpp_error( fatal, 'MPP_GET_PELIST: must first call mpp_domains_init.' )
736  ndivs = size(domain%list(:))
737 
738  if( size(pelist(:)).NE.ndivs ) &
739  call mpp_error( fatal, 'MPP_GET_PELIST: pelist array size does not match domain.' )
740 
741  pelist(:) = domain%list(0:ndivs-1)%pe
742  if( PRESENT(pos) )pos = domain%pos
743  return
744 end subroutine mpp_get_pelist1d
745 
746 !#####################################################################
747 ! <SUBROUTINE NAME="mpp_get_pelist2D" INTERFACE="mpp_get_pelist">
748 ! <IN NAME="domain" TYPE="type(domain2D)"></IN>
749 ! <OUT NAME="pelist" TYPE="integer" DIM="(:)"></OUT>
750 ! <OUT NAME="pos" TYPE="integer"></OUT>
751 ! </SUBROUTINE>
752 subroutine mpp_get_pelist2d( domain, pelist, pos )
753  type(domain2d), intent(in) :: domain
754  integer, intent(out) :: pelist(:)
755  integer, intent(out), optional :: pos
756 
757  if( .NOT.module_is_initialized ) &
758  call mpp_error( fatal, 'MPP_GET_PELIST: must first call mpp_domains_init.' )
759  if( size(pelist(:)).NE.size(domain%list(:)) ) &
760  call mpp_error( fatal, 'MPP_GET_PELIST: pelist array size does not match domain.' )
761 
762  pelist(:) = domain%list(:)%pe
763  if( PRESENT(pos) )pos = domain%pos
764  return
765 end subroutine mpp_get_pelist2d
766 
767 !#####################################################################
768 ! <SUBROUTINE NAME="mpp_get_layout1D" INTERFACE="mpp_get_layout">
769 ! <IN NAME="domain" TYPE="type(domain1D)"></IN>
770 ! <OUT NAME="layout" TYPE="integer"></OUT>
771 ! </SUBROUTINE>
772 subroutine mpp_get_layout1d( domain, layout )
773  type(domain1d), intent(in) :: domain
774  integer, intent(out) :: layout
775 
776  if( .NOT.module_is_initialized ) &
777  call mpp_error( fatal, 'MPP_GET_LAYOUT: must first call mpp_domains_init.' )
778 
779  layout = size(domain%list(:))
780  return
781 end subroutine mpp_get_layout1d
782 
783 !#####################################################################
784 ! <SUBROUTINE NAME="mpp_get_layout2D" INTERFACE="mpp_get_layout">
785 ! <IN NAME="domain" TYPE="type(domain2D)"></IN>
786 ! <OUT NAME="layout" TYPE="integer" DIM="(2)"></OUT>
787 ! </SUBROUTINE>
788 subroutine mpp_get_layout2d( domain, layout )
789  type(domain2d), intent(in) :: domain
790  integer, intent(out) :: layout(2)
791 
792  if( .NOT.module_is_initialized ) &
793  call mpp_error( fatal, 'MPP_GET_LAYOUT: must first call mpp_domains_init.' )
794 
795  layout(1) = size(domain%x(1)%list(:))
796  layout(2) = size(domain%y(1)%list(:))
797  return
798 end subroutine mpp_get_layout2d
799 
800 !#####################################################################
801 
802  !> @brief Returns the shift value in x and y-direction according to domain position..
803  !!
804  !> When domain is symmetry, one extra point maybe needed in
805  !! x- and/or y-direction. This routine will return the shift value based
806  !! on the position
807  !! @code{.F90}
808  !! call mpp_get_domain_shift( domain, ishift, jshift, position )
809  !! @endcode
810 subroutine mpp_get_domain_shift(domain, ishift, jshift, position)
811  type(domain2d), intent(in) :: domain!> predefined data contains 2-d domain decomposition.
812  integer, intent(out) :: ishift, jshift!< return value will be 0 or 1.
813  integer, optional, intent(in) :: position!< position of data. Its value can be CENTER, EAST, NORTH or CORNER.
814  integer :: pos
815 
816  ishift = 0 ; jshift = 0
817  pos = center
818  if(present(position)) pos = position
819 
820  if(domain%symmetry) then ! shift is non-zero only when the domain is symmetry.
821  select case(pos)
822  case(corner)
823  ishift = 1; jshift = 1
824  case(east)
825  ishift = 1
826  case(north)
827  jshift = 1
828  end select
829  end if
830 
831 end subroutine mpp_get_domain_shift
832 
833 !#####################################################################
834 
835  !> Return PE to the righ/left of this PE-domain.
836  subroutine mpp_get_neighbor_pe_1d(domain, direction, pe)
837  type(domain1d), intent(inout) :: domain
838  integer, intent(in) :: direction
839  integer, intent(out) :: pe
840 
841  integer ipos, ipos2, npx
842 
843  pe = null_pe
844  npx = size(domain%list(:)) ! 0..npx-1
845  ipos = domain%pos
846 
847  select case (direction)
848 
849  case (:-1)
850  ! neighbor on the left
851  ipos2 = ipos - 1
852  if(ipos2 < 0) then
853  if(domain%cyclic) then
854  ipos2 = npx-1
855  else
856  ipos2 = -999
857  endif
858  endif
859 
860  case (0)
861  ! identity
862  ipos2 = ipos
863 
864  case (1:)
865  ! neighbor on the right
866  ipos2 = ipos + 1
867  if(ipos2 > npx-1) then
868  if(domain%cyclic) then
869  ipos2 = 0
870  else
871  ipos2 = -999
872  endif
873  endif
874 
875  end select
876 
877  if(ipos2 >= 0) pe = domain%list(ipos2)%pe
878 
879  end subroutine mpp_get_neighbor_pe_1d
880 !#####################################################################
881 
882  !> Return PE North/South/East/West of this PE-domain.
883  !! direction must be NORTH, SOUTH, EAST or WEST.
884  subroutine mpp_get_neighbor_pe_2d(domain, direction, pe)
885  type(domain2d), intent(inout) :: domain
886  integer, intent(in) :: direction
887  integer, intent(out) :: pe
888 
889  integer ipos, jpos, npx, npy, ix, iy, ipos0, jpos0
890 
891  pe = null_pe
892  npx = size(domain%x(1)%list(:)) ! 0..npx-1
893  npy = size(domain%y(1)%list(:)) ! 0..npy-1
894  ipos0 = domain%x(1)%pos
895  jpos0 = domain%y(1)%pos
896 
897  select case (direction)
898  case (north)
899  ix = 0
900  iy = 1
901  case (north_east)
902  ix = 1
903  iy = 1
904  case (east)
905  ix = 1
906  iy = 0
907  case (south_east)
908  ix = 1
909  iy =-1
910  case (south)
911  ix = 0
912  iy =-1
913  case (south_west)
914  ix =-1
915  iy =-1
916  case (west)
917  ix =-1
918  iy = 0
919  case (north_west)
920  ix =-1
921  iy = 1
922 
923  case default
924  call mpp_error( fatal, &
925  & 'MPP_GET_NEIGHBOR_PE_2D: direction must be either NORTH, ' &
926  & // 'SOUTH, EAST, WEST, NORTH_EAST, SOUTH_EAST, SOUTH_WEST or NORTH_WEST')
927  end select
928 
929  ipos = ipos0 + ix
930  jpos = jpos0 + iy
931 
932 
933  if( (ipos < 0 .or. ipos > npx-1) .and. domain%x(1)%cyclic ) then
934  ! E/W cyclic domain
935  ipos = modulo(ipos, npx)
936  endif
937 
938  if( (ipos < 0 .and. btest(domain%fold,west)) .or. &
939  & (ipos > npx-1 .and. btest(domain%fold,east)) ) then
940  ! E or W folded domain
941  ipos = ipos0
942  jpos = npy-jpos-1
943  endif
944 
945  if( (jpos < 0 .or. jpos > npy-1) .and. domain%y(1)%cyclic ) then
946  ! N/S cyclic
947  jpos = modulo(jpos, npy)
948  endif
949 
950  if( (jpos < 0 .and. btest(domain%fold,south)) .or. &
951  & (jpos > npy-1 .and. btest(domain%fold,north)) ) then
952  ! N or S folded
953  ipos = npx-ipos-1
954  jpos = jpos0
955  endif
956 
957  ! get the PE number
958  pe = null_pe
959  if(ipos >= 0 .and. ipos <= npx-1 .and. jpos >= 0 .and. jpos <= npy-1) then
960  pe = domain%pearray(ipos, jpos)
961  endif
962 
963 
964  end subroutine mpp_get_neighbor_pe_2d
965 
966 
967 !#######################################################################
968 
969  subroutine nullify_domain2d_list(domain)
970  type(domain2d), intent(inout) :: domain
971 
972  domain%list =>null()
973 
974  end subroutine nullify_domain2d_list
975 
976 !#######################################################################
977  function mpp_domain_is_symmetry(domain)
978  type(domain2d), intent(in) :: domain
979  logical :: mpp_domain_is_symmetry
980 
981  mpp_domain_is_symmetry = domain%symmetry
982  return
983 
984  end function mpp_domain_is_symmetry
985 
986 !#######################################################################
987  function mpp_domain_is_initialized(domain)
988  type(domain2d), intent(in) :: domain
989  logical :: mpp_domain_is_initialized
990 
991  mpp_domain_is_initialized = domain%initialized
992 
993  return
994 
995  end function mpp_domain_is_initialized
996 
997 !#######################################################################
998  !--- private routine used only for mpp_update_domains. This routine will
999  !--- compare whalo, ehalo, shalo, nhalo with the halo size when defining "domain"
1000  !--- to decide if update is needed. Also it check the sign of whalo, ehalo, shalo and nhalo.
1001  function domain_update_is_needed(domain, whalo, ehalo, shalo, nhalo)
1002  type(domain2d), intent(in) :: domain
1003  integer, intent(in) :: whalo, ehalo, shalo, nhalo
1004  logical :: domain_update_is_needed
1005 
1006  domain_update_is_needed = .true.
1007 
1008  if(whalo == 0 .AND. ehalo==0 .AND. shalo == 0 .AND. nhalo==0 ) then
1009  domain_update_is_needed = .false.
1010  if( debug )call mpp_error(note, &
1011  'mpp_domains_util.inc: halo size to be updated are all zero, no update will be done')
1012  return
1013  end if
1014  if( (whalo == -domain%whalo .AND. domain%whalo .NE. 0) .or. &
1015  (ehalo == -domain%ehalo .AND. domain%ehalo .NE. 0) .or. &
1016  (shalo == -domain%shalo .AND. domain%shalo .NE. 0) .or. &
1017  (nhalo == -domain%nhalo .AND. domain%nhalo .NE. 0) ) then
1018  domain_update_is_needed = .false.
1019  call mpp_error(note, 'mpp_domains_util.inc: at least one of w/e/s/n halo size to be updated '// &
1020  'is the inverse of the original halo when defining domain, no update will be done')
1021  return
1022  end if
1023 
1024  end function domain_update_is_needed
1025 !#######################################################################
1026  !> this routine found the domain has the same halo size with the input
1027  !! whalo, ehalo,
1028  function search_update_overlap(domain, whalo, ehalo, shalo, nhalo, position)
1029  type(domain2d), intent(inout) :: domain
1030  integer, intent(in) :: whalo, ehalo, shalo, nhalo
1031  integer, intent(in) :: position
1032  type(overlapspec), pointer :: search_update_overlap
1033  type(overlapspec), pointer :: update_ref
1034  type(overlapspec), pointer :: check => null()
1035  integer :: ishift, jshift, shift
1036 
1037  shift = 0; if(domain%symmetry) shift = 1
1038  select case(position)
1039  case (center)
1040  update_ref => domain%update_T
1041  ishift = 0; jshift = 0
1042  case (corner)
1043  update_ref => domain%update_C
1044  ishift = shift; jshift = shift
1045  case (north)
1046  update_ref => domain%update_N
1047  ishift = 0; jshift = shift
1048  case (east)
1049  update_ref => domain%update_E
1050  ishift = shift; jshift = 0
1051  case default
1052  call mpp_error(fatal,"mpp_domains_util.inc(search_update_overlap): position should be CENTER|CORNER|EAST|NORTH")
1053  end select
1054 
1055  search_update_overlap => update_ref
1056 
1057  do
1058  if(whalo == search_update_overlap%whalo .AND. ehalo == search_update_overlap%ehalo .AND. &
1059  shalo == search_update_overlap%shalo .AND. nhalo == search_update_overlap%nhalo ) then
1060  exit ! found domain
1061  endif
1062  !--- if not found, switch to next
1063  if(.NOT. ASSOCIATED(search_update_overlap%next)) then
1064  allocate(search_update_overlap%next)
1066  if(domain%fold .NE. 0) then
1067  call compute_overlaps(domain, position, search_update_overlap, check, &
1068  ishift, jshift, 0, 0, whalo, ehalo, shalo, nhalo)
1069  else
1070  call set_overlaps(domain, update_ref, search_update_overlap, whalo, ehalo, shalo, nhalo )
1071  endif
1072  exit
1073  else
1075  end if
1076 
1077  end do
1078 
1079  update_ref => null()
1080 
1081  end function search_update_overlap
1082 
1083 !#######################################################################
1084  !> this routine finds the check at certain position
1085  function search_check_overlap(domain, position)
1086  type(domain2d), intent(in) :: domain
1087  integer, intent(in) :: position
1088  type(overlapspec), pointer :: search_check_overlap
1089 
1090  select case(position)
1091  case (center)
1092  search_check_overlap => null()
1093  case (corner)
1094  search_check_overlap => domain%check_C
1095  case (north)
1096  search_check_overlap => domain%check_N
1097  case (east)
1098  search_check_overlap => domain%check_E
1099  case default
1100  call mpp_error(fatal,"mpp_domains_util.inc(search_check_overlap): position should be CENTER|CORNER|EAST|NORTH")
1101  end select
1102 
1103  end function search_check_overlap
1104 
1105 !#######################################################################
1106  !> This routine finds the bound at certain position
1107  function search_bound_overlap(domain, position)
1108  type(domain2d), intent(in) :: domain
1109  integer, intent(in) :: position
1110  type(overlapspec), pointer :: search_bound_overlap
1111 
1112  select case(position)
1113  case (center)
1114  search_bound_overlap => null()
1115  case (corner)
1116  search_bound_overlap => domain%bound_C
1117  case (north)
1118  search_bound_overlap => domain%bound_N
1119  case (east)
1120  search_bound_overlap => domain%bound_E
1121  case default
1122  call mpp_error(fatal,"mpp_domains_util.inc(search_bound_overlap): position should be CENTER|CORNER|EAST|NORTH")
1123  end select
1124 
1125  end function search_bound_overlap
1126 
1127  !########################################################################
1128  !> Returns the tile_id on current pe
1129  function mpp_get_tile_id(domain)
1130  type(domain2d), intent(in) :: domain
1131  integer, dimension(size(domain%tile_id(:))) :: mpp_get_tile_id
1132 
1133  mpp_get_tile_id = domain%tile_id
1134  return
1135 
1136  end function mpp_get_tile_id
1137 
1138  !#######################################################################
1139  !> Return the tile_id on current pelist. one-tile-per-pe is assumed.
1140  subroutine mpp_get_tile_list(domain, tiles)
1141  type(domain2d), intent(in) :: domain
1142  integer, intent(inout) :: tiles(:)
1143  integer :: i
1144 
1145  if( size(tiles(:)).NE.size(domain%list(:)) ) &
1146  call mpp_error( fatal, 'mpp_get_tile_list: tiles array size does not match domain.' )
1147  do i = 1, size(tiles(:))
1148  if(size(domain%list(i-1)%tile_id(:)) > 1) call mpp_error( fatal, &
1149  'mpp_get_tile_list: only support one-tile-per-pe now, contact developer');
1150  tiles(i) = domain%list(i-1)%tile_id(1)
1151  end do
1152 
1153  end subroutine mpp_get_tile_list
1154 
1155  !########################################################################
1156  !> Returns number of tiles in mosaic
1157  function mpp_get_ntile_count(domain)
1158  type(domain2d), intent(in) :: domain
1159  integer :: mpp_get_ntile_count
1160 
1161  mpp_get_ntile_count = domain%ntiles
1162  return
1163 
1164  end function mpp_get_ntile_count
1165 
1166  !########################################################################
1167  !> Returns number of tile on current pe
1168  function mpp_get_current_ntile(domain)
1169  type(domain2d), intent(in) :: domain
1170  integer :: mpp_get_current_ntile
1171 
1172  mpp_get_current_ntile = size(domain%tile_id(:))
1173  return
1174 
1175  end function mpp_get_current_ntile
1176 
1177  !#######################################################################
1178  !> Returns if current pe is the root pe of the tile, if number of tiles on current pe
1179  !! is greater than 1, will return true, if isc==isg and jsc==jsg also will return true,
1180  !! otherwise false will be returned.
1182  type(domain2d), intent(in) :: domain
1183  logical :: mpp_domain_is_tile_root_pe
1184 
1185  mpp_domain_is_tile_root_pe = domain%pe == domain%tile_root_pe;
1186 
1187  end function mpp_domain_is_tile_root_pe
1188 
1189  !#########################################################################
1190  !> Returns number of processors used on current tile.
1191  function mpp_get_tile_npes(domain)
1192  type(domain2d), intent(in) :: domain
1193  integer :: mpp_get_tile_npes
1194  integer :: i, tile
1195 
1196  !--- When there is more than one tile on this pe, we assume each tile will be
1197  !--- limited to this pe.
1198  if(size(domain%tile_id(:)) > 1) then
1199  mpp_get_tile_npes = 1
1200  else
1201  mpp_get_tile_npes = 0
1202  tile = domain%tile_id(1)
1203  do i = 0, size(domain%list(:))-1
1204  if(tile == domain%list(i)%tile_id(1) ) mpp_get_tile_npes = mpp_get_tile_npes + 1
1205  end do
1206  endif
1207 
1208  end function mpp_get_tile_npes
1209 
1210  !########################################################################
1211  !> Get the processors list used on current tile.
1212  subroutine mpp_get_tile_pelist(domain, pelist)
1213  type(domain2d), intent(in) :: domain
1214  integer, intent(inout) :: pelist(:)
1215  integer :: npes_on_tile
1216  integer :: i, tile, pos
1217 
1218  npes_on_tile = mpp_get_tile_npes(domain)
1219  if(size(pelist(:)) .NE. npes_on_tile) call mpp_error(fatal, &
1220  "mpp_domains_util.inc(mpp_get_tile_pelist): size(pelist) does not equal npes on current tile")
1221  tile = domain%tile_id(1)
1222  pos = 0
1223  do i = 0, size(domain%list(:))-1
1224  if(tile == domain%list(i)%tile_id(1)) then
1225  pos = pos+1
1226  pelist(pos) = domain%list(i)%pe
1227  endif
1228  enddo
1229 
1230  return
1231 
1232  end subroutine mpp_get_tile_pelist
1233 
1234 !#####################################################################
1235 subroutine mpp_get_tile_compute_domains( domain, xbegin, xend, ybegin, yend, position )
1236  type(domain2d), intent(in) :: domain
1237  integer, intent(out), dimension(:) :: xbegin, xend, ybegin, yend
1238  integer, intent(in ), optional :: position
1239 
1240  integer :: i, ishift, jshift
1241  integer :: npes_on_tile, pos, tile
1242 
1243  call mpp_get_domain_shift( domain, ishift, jshift, position )
1244 
1245 
1246  if( .NOT.module_is_initialized ) &
1247  call mpp_error( fatal, 'mpp_get_compute_domains2D: must first call mpp_domains_init.' )
1248 
1249  npes_on_tile = mpp_get_tile_npes(domain)
1250  if(size(xbegin(:)) .NE. npes_on_tile) call mpp_error(fatal, &
1251  "mpp_domains_util.inc(mpp_get_compute_domains2D): size(xbegin) does not equal npes on current tile")
1252  if(size(xend(:)) .NE. npes_on_tile) call mpp_error(fatal, &
1253  "mpp_domains_util.inc(mpp_get_compute_domains2D): size(xend) does not equal npes on current tile")
1254  if(size(ybegin(:)) .NE. npes_on_tile) call mpp_error(fatal, &
1255  "mpp_domains_util.inc(mpp_get_compute_domains2D): size(ybegin) does not equal npes on current tile")
1256  if(size(yend(:)) .NE. npes_on_tile) call mpp_error(fatal, &
1257  "mpp_domains_util.inc(mpp_get_compute_domains2D): size(yend) does not equal npes on current tile")
1258 
1259  tile = domain%tile_id(1)
1260  pos = 0
1261  do i = 0, size(domain%list(:))-1
1262  if(tile == domain%list(i)%tile_id(1)) then
1263  pos = pos+1
1264  xbegin(pos) = domain%list(i)%x(1)%compute%begin
1265  xend(pos) = domain%list(i)%x(1)%compute%end + ishift
1266  ybegin(pos) = domain%list(i)%y(1)%compute%begin
1267  yend(pos) = domain%list(i)%y(1)%compute%end + jshift
1268  endif
1269  enddo
1270 
1271  return
1272 
1273 end subroutine mpp_get_tile_compute_domains
1274 
1275 
1276 
1277  !#############################################################################
1278  function mpp_get_num_overlap(domain, action, p, position)
1279  type(domain2d), intent(in) :: domain
1280  integer, intent(in) :: action
1281  integer, intent(in) :: p
1282  integer, optional, intent(in) :: position
1283  integer :: mpp_get_num_overlap
1284  type(overlapspec), pointer :: update => null()
1285  integer :: pos
1286 
1287  pos = center
1288  if(present(position)) pos = position
1289  select case(pos)
1290  case (center)
1291  update => domain%update_T
1292  case (corner)
1293  update => domain%update_C
1294  case (east)
1295  update => domain%update_E
1296  case (north)
1297  update => domain%update_N
1298  case default
1299  call mpp_error( fatal, "mpp_domains_mod(mpp_get_num_overlap): invalid option of position")
1300  end select
1301 
1302  if(action == event_send) then
1303  if(p< 1 .OR. p > update%nsend) call mpp_error( fatal, &
1304  "mpp_domains_mod(mpp_get_num_overlap): p should be between 1 and update%nsend")
1305  mpp_get_num_overlap = update%send(p)%count
1306  else if(action == event_recv) then
1307  if(p< 1 .OR. p > update%nrecv) call mpp_error( fatal, &
1308  "mpp_domains_mod(mpp_get_num_overlap): p should be between 1 and update%nrecv")
1309  mpp_get_num_overlap = update%recv(p)%count
1310  else
1311  call mpp_error( fatal, "mpp_domains_mod(mpp_get_num_overlap): invalid option of action")
1312  end if
1313 
1314  end function mpp_get_num_overlap
1315 
1316  !#############################################################################
1317  subroutine mpp_get_update_size(domain, nsend, nrecv, position)
1318  type(domain2d), intent(in) :: domain
1319  integer, intent(out) :: nsend, nrecv
1320  integer, optional, intent(in) :: position
1321  integer :: pos
1322 
1323  pos = center
1324  if(present(position)) pos = position
1325  select case(pos)
1326  case (center)
1327  nsend = domain%update_T%nsend
1328  nrecv = domain%update_T%nrecv
1329  case (corner)
1330  nsend = domain%update_C%nsend
1331  nrecv = domain%update_C%nrecv
1332  case (east)
1333  nsend = domain%update_E%nsend
1334  nrecv = domain%update_E%nrecv
1335  case (north)
1336  nsend = domain%update_N%nsend
1337  nrecv = domain%update_N%nrecv
1338  case default
1339  call mpp_error( fatal, "mpp_domains_mod(mpp_get_update_size): invalid option of position")
1340  end select
1341 
1342  end subroutine mpp_get_update_size
1343 
1344  !#############################################################################
1345  subroutine mpp_get_update_pelist(domain, action, pelist, position)
1346  type(domain2d), intent(in) :: domain
1347  integer, intent(in) :: action
1348  integer, intent(inout) :: pelist(:)
1349  integer, optional, intent(in) :: position
1350  type(overlapspec), pointer :: update => null()
1351  integer :: pos, p
1352 
1353  pos = center
1354  if(present(position)) pos = position
1355  select case(pos)
1356  case (center)
1357  update => domain%update_T
1358  case (corner)
1359  update => domain%update_C
1360  case (east)
1361  update => domain%update_E
1362  case (north)
1363  update => domain%update_N
1364  case default
1365  call mpp_error( fatal, "mpp_domains_mod(mpp_get_update_pelist): invalid option of position")
1366  end select
1367 
1368  if(action == event_send) then
1369  if(size(pelist) .NE. update%nsend) call mpp_error( fatal, &
1370  "mpp_domains_mod(mpp_get_update_pelist): size of pelist does not match update%nsend")
1371  do p = 1, update%nsend
1372  pelist(p) = update%send(p)%pe
1373  enddo
1374  else if(action == event_recv) then
1375  if(size(pelist) .NE. update%nrecv) call mpp_error( fatal, &
1376  "mpp_domains_mod(mpp_get_update_pelist): size of pelist does not match update%nrecv")
1377  do p = 1, update%nrecv
1378  pelist(p) = update%recv(p)%pe
1379  enddo
1380  else
1381  call mpp_error( fatal, "mpp_domains_mod(mpp_get_update_pelist): invalid option of action")
1382  end if
1383 
1384  end subroutine mpp_get_update_pelist
1385 
1386  !#############################################################################
1387  subroutine mpp_get_overlap(domain, action, p, is, ie, js, je, dir, rot, position)
1388  type(domain2d), intent(in) :: domain
1389  integer, intent(in) :: action
1390  integer, intent(in) :: p
1391  integer, dimension(:), intent(out) :: is, ie, js, je
1392  integer, dimension(:), intent(out) :: dir, rot
1393  integer, optional, intent(in) :: position
1394  type(overlapspec), pointer :: update => null()
1395  type(overlap_type), pointer :: overlap => null()
1396  integer :: count, pos
1397 
1398  pos = center
1399  if(present(position)) pos = position
1400  select case(pos)
1401  case (center)
1402  update => domain%update_T
1403  case (corner)
1404  update => domain%update_C
1405  case (east)
1406  update => domain%update_E
1407  case (north)
1408  update => domain%update_N
1409  case default
1410  call mpp_error( fatal, "mpp_domains_mod(mpp_get_overlap): invalid option of position")
1411  end select
1412 
1413  if(action == event_send) then
1414  overlap => update%send(p)
1415  else if(action == event_recv) then
1416  overlap => update%recv(p)
1417  else
1418  call mpp_error( fatal, "mpp_domains_mod(mpp_get_overlap): invalid option of action")
1419  end if
1420 
1421  count = overlap%count
1422  if(size(is(:)) .NE. count .OR. size(ie(:)) .NE. count .OR. size(js(:)) .NE. count .OR. &
1423  size(je(:)) .NE. count .OR. size(dir(:)) .NE. count .OR. size(rot(:)) .NE. count ) &
1424  call mpp_error( fatal, &
1425  & "mpp_domains_mod(mpp_get_overlap): size mismatch between number of overlap and array size")
1426 
1427  is = overlap%is (1:count)
1428  ie = overlap%ie (1:count)
1429  js = overlap%js (1:count)
1430  je = overlap%je (1:count)
1431  dir = overlap%dir (1:count)
1432  rot = overlap%rotation(1:count)
1433 
1434  update => null()
1435  overlap => null()
1436 
1437  end subroutine mpp_get_overlap
1438 
1439  !##################################################################
1440  function mpp_get_domain_name(domain)
1441  type(domain2d), intent(in) :: domain
1442  character(len=NAME_LENGTH) :: mpp_get_domain_name
1443 
1444  mpp_get_domain_name = domain%name
1445 
1446  end function mpp_get_domain_name
1447 
1448  !#################################################################
1449  function mpp_get_domain_root_pe(domain)
1450  type(domain2d), intent(in) :: domain
1451  integer :: mpp_get_domain_root_pe
1452 
1453  mpp_get_domain_root_pe = domain%list(0)%pe
1454 
1455  end function mpp_get_domain_root_pe
1456  !#################################################################
1457  function mpp_get_domain_npes(domain)
1458  type(domain2d), intent(in) :: domain
1459  integer :: mpp_get_domain_npes
1460 
1461  mpp_get_domain_npes = size(domain%list(:))
1462 
1463  return
1464 
1465  end function mpp_get_domain_npes
1466 
1467  !################################################################
1468  subroutine mpp_get_domain_pelist(domain, pelist)
1469  type(domain2d), intent(in) :: domain
1470  integer, intent(out) :: pelist(:)
1471  integer :: p
1472 
1473  if(size(pelist(:)) .NE. size(domain%list(:)) ) then
1474  call mpp_error(fatal, .NE."mpp_get_domain_pelist: size(pelist(:)) size(domain%list(:)) ")
1475  endif
1476 
1477  do p = 0, size(domain%list(:))-1
1478  pelist(p+1) = domain%list(p)%pe
1479  enddo
1480 
1481  return
1482 
1483  end subroutine mpp_get_domain_pelist
1484 
1485  !#################################################################
1486  function mpp_get_io_domain_layout(domain)
1487  type(domain2d), intent(in) :: domain
1488  integer, dimension(2) :: mpp_get_io_domain_layout
1489 
1490  mpp_get_io_domain_layout = domain%io_layout
1491 
1492  end function mpp_get_io_domain_layout
1493 
1494  !################################################################
1495  function get_rank_send(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
1496  type(domain2d), intent(in) :: domain
1497  type(overlapspec), intent(in) :: overlap_x, overlap_y
1498  integer, intent(out) :: rank_x, rank_y, ind_x, ind_y
1499  integer :: get_rank_send
1500  integer :: nlist, nsend_x, nsend_y
1501 
1502  nlist = size(domain%list(:))
1503  nsend_x = overlap_x%nsend
1504  nsend_y = overlap_y%nsend
1505  rank_x = nlist+1
1506  rank_y = nlist+1
1507  if(nsend_x>0) rank_x = overlap_x%send(1)%pe - domain%pe
1508  if(nsend_y>0) rank_y = overlap_y%send(1)%pe - domain%pe
1509  if(rank_x .LT. 0) rank_x = rank_x + nlist
1510  if(rank_y .LT. 0) rank_y = rank_y + nlist
1511  get_rank_send = min(rank_x, rank_y)
1512  ind_x = nsend_x + 1
1513  ind_y = nsend_y + 1
1514  if(get_rank_send < nlist+1) then
1515  if(nsend_x>0) ind_x = 1
1516  if(nsend_y>0) ind_y = 1
1517  endif
1518 
1519  end function get_rank_send
1520 
1521  !############################################################################
1522  function get_rank_recv(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
1523  type(domain2d), intent(in) :: domain
1524  type(overlapspec), intent(in) :: overlap_x, overlap_y
1525  integer, intent(out) :: rank_x, rank_y, ind_x, ind_y
1526  integer :: get_rank_recv
1527  integer :: nlist, nrecv_x, nrecv_y
1528 
1529  nlist = size(domain%list(:))
1530  nrecv_x = overlap_x%nrecv
1531  nrecv_y = overlap_y%nrecv
1532  rank_x = -1
1533  rank_y = -1
1534  if(nrecv_x>0) then
1535  rank_x = overlap_x%recv(1)%pe - domain%pe
1536  if(rank_x .LE. 0) rank_x = rank_x + nlist
1537  endif
1538  if(nrecv_y>0) then
1539  rank_y = overlap_y%recv(1)%pe - domain%pe
1540  if(rank_y .LE. 0) rank_y = rank_y + nlist
1541  endif
1542  get_rank_recv = max(rank_x, rank_y)
1543  ind_x = nrecv_x + 1
1544  ind_y = nrecv_y + 1
1545  if(get_rank_recv < nlist+1) then
1546  if(nrecv_x>0) ind_x = 1
1547  if(nrecv_y>0) ind_y = 1
1548  endif
1549 
1550  end function get_rank_recv
1551 
1552  function get_vector_recv(domain, update_x, update_y, ind_x, ind_y, start_pos, pelist)
1553  type(domain2d), intent(in) :: domain
1554  type(overlapspec), intent(in) :: update_x, update_y
1555  integer, intent(out) :: ind_x(:), ind_y(:)
1556  integer, intent(out) :: start_pos(:)
1557  integer, intent(out) :: pelist(:)
1558  integer :: nlist, nrecv_x, nrecv_y, ntot, n
1559  integer :: ix, iy, rank_x, rank_y, cur_pos
1560  integer :: get_vector_recv
1561 
1562  nlist = size(domain%list(:))
1563  nrecv_x = update_x%nrecv
1564  nrecv_y = update_y%nrecv
1565 
1566  ntot = nrecv_x + nrecv_y
1567 
1568  n = 1
1569  ix = 1
1570  iy = 1
1571  ind_x = -1
1572  ind_y = -1
1573  get_vector_recv = 0
1574  cur_pos = 0
1575  do while (n<=ntot)
1576  if(ix <= nrecv_x ) then
1577  rank_x = update_x%recv(ix)%pe-domain%pe
1578  if(rank_x .LE. 0) rank_x = rank_x + nlist
1579  else
1580  rank_x = -1
1581  endif
1582  if(iy <= nrecv_y ) then
1583  rank_y = update_y%recv(iy)%pe-domain%pe
1584  if(rank_y .LE. 0) rank_y = rank_y + nlist
1585  else
1586  rank_y = -1
1587  endif
1589  start_pos(get_vector_recv) = cur_pos
1590  if( rank_x == rank_y ) then
1591  n = n+2
1592  ind_x(get_vector_recv) = ix
1593  ind_y(get_vector_recv) = iy
1594  cur_pos = cur_pos + update_x%recv(ix)%totsize + update_y%recv(iy)%totsize
1595  pelist(get_vector_recv) = update_x%recv(ix)%pe
1596  ix = ix + 1
1597  iy = iy + 1
1598  else if ( rank_x > rank_y ) then
1599  n = n+1
1600  ind_x(get_vector_recv) = ix
1601  ind_y(get_vector_recv) = -1
1602  cur_pos = cur_pos + update_x%recv(ix)%totsize
1603  pelist(get_vector_recv) = update_x%recv(ix)%pe
1604  ix = ix + 1
1605  else if ( rank_y > rank_x ) then
1606  n = n+1
1607  ind_x(get_vector_recv) = -1
1608  ind_y(get_vector_recv) = iy
1609  cur_pos = cur_pos + update_y%recv(iy)%totsize
1610  pelist(get_vector_recv) = update_y%recv(iy)%pe
1611  iy = iy+1
1612  endif
1613  end do
1614 
1615 
1616  end function get_vector_recv
1617 
1618  function get_vector_send(domain, update_x, update_y, ind_x, ind_y, start_pos, pelist)
1619  type(domain2d), intent(in) :: domain
1620  type(overlapspec), intent(in) :: update_x, update_y
1621  integer, intent(out) :: ind_x(:), ind_y(:)
1622  integer, intent(out) :: start_pos(:)
1623  integer, intent(out) :: pelist(:)
1624  integer :: nlist, nsend_x, nsend_y, ntot, n
1625  integer :: ix, iy, rank_x, rank_y, cur_pos
1626  integer :: get_vector_send
1627 
1628  nlist = size(domain%list(:))
1629  nsend_x = update_x%nsend
1630  nsend_y = update_y%nsend
1631 
1632  ntot = nsend_x + nsend_y
1633  n = 1
1634  ix = 1
1635  iy = 1
1636  ind_x = -1
1637  ind_y = -1
1638  get_vector_send = 0
1639  cur_pos = 0
1640  do while (n<=ntot)
1641  if(ix <= nsend_x ) then
1642  rank_x = update_x%send(ix)%pe-domain%pe
1643  if(rank_x .LT. 0) rank_x = rank_x + nlist
1644  else
1645  rank_x = nlist+1
1646  endif
1647  if(iy <= nsend_y ) then
1648  rank_y = update_y%send(iy)%pe-domain%pe
1649  if(rank_y .LT. 0) rank_y = rank_y + nlist
1650  else
1651  rank_y = nlist+1
1652  endif
1654  start_pos(get_vector_send) = cur_pos
1655 
1656  if( rank_x == rank_y ) then
1657  n = n+2
1658  ind_x(get_vector_send) = ix
1659  ind_y(get_vector_send) = iy
1660  cur_pos = cur_pos + update_x%send(ix)%totsize + update_y%send(iy)%totsize
1661  pelist(get_vector_send) = update_x%send(ix)%pe
1662  ix = ix + 1
1663  iy = iy + 1
1664  else if ( rank_x < rank_y ) then
1665  n = n+1
1666  ind_x(get_vector_send) = ix
1667  ind_y(get_vector_send) = -1
1668  cur_pos = cur_pos + update_x%send(ix)%totsize
1669  pelist(get_vector_send) = update_x%send(ix)%pe
1670  ix = ix + 1
1671  else if ( rank_y < rank_x ) then
1672  n = n+1
1673  ind_x(get_vector_send) = -1
1674  ind_y(get_vector_send) = iy
1675  cur_pos = cur_pos + update_y%send(iy)%totsize
1676  pelist(get_vector_send) = update_y%send(iy)%pe
1677  iy = iy+1
1678  endif
1679  end do
1680 
1681 
1682  end function get_vector_send
1683 
1684 
1685  !############################################################################
1686  function get_rank_unpack(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
1687  type(domain2d), intent(in) :: domain
1688  type(overlapspec), intent(in) :: overlap_x, overlap_y
1689  integer, intent(out) :: rank_x, rank_y, ind_x, ind_y
1690  integer :: get_rank_unpack
1691  integer :: nlist, nrecv_x, nrecv_y
1692 
1693  nlist = size(domain%list(:))
1694  nrecv_x = overlap_x%nrecv
1695  nrecv_y = overlap_y%nrecv
1696 
1697  rank_x = nlist+1
1698  rank_y = nlist+1
1699  if(nrecv_x>0) rank_x = overlap_x%recv(nrecv_x)%pe - domain%pe
1700  if(nrecv_y>0) rank_y = overlap_y%recv(nrecv_y)%pe - domain%pe
1701  if(rank_x .LE.0) rank_x = rank_x + nlist
1702  if(rank_y .LE.0) rank_y = rank_y + nlist
1703 
1704  get_rank_unpack = min(rank_x, rank_y)
1705  ind_x = 0
1706  ind_y = 0
1707  if(get_rank_unpack < nlist+1) then
1708  if(nrecv_x >0) ind_x = nrecv_x
1709  if(nrecv_y >0) ind_y = nrecv_y
1710  endif
1711 
1712  end function get_rank_unpack
1713 
1714  function get_mesgsize(overlap, do_dir)
1715  type(overlap_type), intent(in) :: overlap
1716  logical, intent(in) :: do_dir(:)
1717  integer :: get_mesgsize
1718  integer :: n, dir
1719 
1720  get_mesgsize = 0
1721  do n = 1, overlap%count
1722  dir = overlap%dir(n)
1723  if(do_dir(dir)) then
1724  get_mesgsize = get_mesgsize + overlap%msgsize(n)
1725  end if
1726  end do
1727 
1728  end function get_mesgsize
1729 
1730  !#############################################################################
1731  subroutine mpp_set_domain_symmetry(domain, symmetry)
1732  type(domain2d), intent(inout) :: domain
1733  logical, intent(in ) :: symmetry
1734 
1735  domain%symmetry = symmetry
1736 
1737  end subroutine mpp_set_domain_symmetry
1738 
1739  !> @brief Copies input 1d domain to the output 1d domain
1740  recursive subroutine mpp_copy_domain1d(domain_in, domain_out)
1741  type(domain1d), intent(in) :: domain_in !< Input domain
1742  type(domain1d), intent(inout) :: domain_out !< Output domain
1743 
1744  integer :: i !< For loop
1745  integer :: starting !< Starting bounds
1746  integer :: ending !< Ending bounds
1747 
1748  domain_out%compute = domain_in%compute
1749  domain_out%domain_data = domain_in%domain_data
1750  domain_out%global = domain_in%global
1751  domain_out%memory = domain_in%memory
1752  domain_out%cyclic = domain_in%cyclic
1753  domain_out%pe = domain_in%pe
1754  domain_out%pos = domain_in%pos
1755 
1756  if (associated(domain_in%list)) then
1757  starting = lbound(domain_in%list, 1)
1758  ending = ubound(domain_in%list, 1)
1759  if (associated(domain_out%list)) deallocate(domain_out%list) !< Check if allocated
1760  allocate(domain_out%list(starting:ending))
1761 
1762  do i = starting, ending
1763  call mpp_copy_domain1d(domain_in%list(i), domain_out%list(i))
1764  enddo
1765 
1766  endif
1767 
1768  end subroutine mpp_copy_domain1d
1769 
1770  !#################################################################
1771  !> @brief Copies input 2d domain to the output 2d domain
1772  subroutine mpp_copy_domain2d(domain_in, domain_out)
1773  type(domain2d), intent(in) :: domain_in !< Input domain
1774  type(domain2d), intent(inout) :: domain_out !< Output domain
1775 
1776  integer :: n,i !< For loops
1777  integer :: ntiles !< Number of tiles
1778  integer :: starting(2) !< Starting bounds
1779  integer :: ending(2) !< Ending bounds
1780 
1781  if (associated(domain_out%x)) then
1782  call mpp_error(fatal, "mpp_copy_domain: domain_out is already set")
1783  endif
1784 
1785  domain_out%id = domain_in%id
1786  domain_out%pe = domain_in%pe
1787  domain_out%fold = domain_in%fold
1788  domain_out%pos = domain_in%pos
1789  domain_out%symmetry = domain_in%symmetry
1790  domain_out%whalo = domain_in%whalo
1791  domain_out%ehalo = domain_in%ehalo
1792  domain_out%shalo = domain_in%shalo
1793  domain_out%nhalo = domain_in%nhalo
1794  domain_out%ntiles = domain_in%ntiles
1795  domain_out%max_ntile_pe = domain_in%max_ntile_pe
1796  domain_out%ncontacts = domain_in%ncontacts
1797  domain_out%rotated_ninety = domain_in%rotated_ninety
1798  domain_out%initialized = domain_in%initialized
1799  domain_out%tile_root_pe = domain_in%tile_root_pe
1800  domain_out%io_layout = domain_in%io_layout
1801  domain_out%name = domain_in%name
1802 
1803  ntiles = size(domain_in%x(:))
1804  allocate(domain_out%x(ntiles), domain_out%y(ntiles), domain_out%tile_id(ntiles) )
1805  do n = 1, ntiles
1806  call mpp_copy_domain1d(domain_in%x(n), domain_out%x(n))
1807  call mpp_copy_domain1d(domain_in%y(n), domain_out%y(n))
1808  enddo
1809 
1810  if (associated(domain_in%pearray)) then
1811  starting = lbound(domain_in%pearray)
1812  ending = ubound(domain_in%pearray)
1813 
1814  allocate(domain_out%pearray(starting(1):ending(1), starting(2):ending(2)))
1815  domain_out%pearray=domain_in%pearray
1816  endif
1817 
1818  if (associated(domain_in%tile_id)) then
1819  starting(1) = lbound(domain_in%tile_id,1)
1820  ending(1) = ubound(domain_in%tile_id,1)
1821 
1822  allocate(domain_out%tile_id(starting(1):ending(1)))
1823  domain_out%tile_id = domain_in%tile_id
1824  endif
1825 
1826  if (associated(domain_in%tile_id_all)) then
1827  starting(1) = lbound(domain_in%tile_id_all,1)
1828  ending(1) = ubound(domain_in%tile_id_all,1)
1829 
1830  allocate(domain_out%tile_id_all(starting(1):ending(1)))
1831  domain_out%tile_id_all = domain_in%tile_id_all
1832  endif
1833 
1834  if (associated(domain_in%list)) then
1835  starting(1) = lbound(domain_in%list,1)
1836  ending(1) = ubound(domain_in%list,1)
1837 
1838  allocate(domain_out%list(starting(1):ending(1)))
1839  do i = starting(1), ending(1)
1840  call mpp_copy_domain2d_spec(domain_in%list(i),domain_out%list(i))
1841  enddo
1842  endif
1843 
1844  return
1845 
1846  end subroutine mpp_copy_domain2d
1847 
1848  !> @brief Copies input 2d domain spec to the output 2d domain spec
1849  subroutine mpp_copy_domain2d_spec(domain2D_spec_in, domain2d_spec_out)
1850  type(domain2d_spec), intent(in) :: domain2D_spec_in !< Input
1851  type(domain2d_spec), intent(out) :: domain2D_spec_out !< Output
1852 
1853  integer :: starting !< Starting bounds
1854  integer :: ending !< Ending bounds
1855  integer :: i !< For loop
1856 
1857  domain2d_spec_out%pe = domain2d_spec_in%pe
1858  domain2d_spec_out%pos = domain2d_spec_in%pos
1859  domain2d_spec_out%tile_root_pe = domain2d_spec_in%tile_root_pe
1860 
1861  if (associated(domain2d_spec_in%tile_id)) then
1862  starting = lbound(domain2d_spec_in%tile_id,1)
1863  ending = ubound(domain2d_spec_in%tile_id,1)
1864 
1865  if (associated(domain2d_spec_out%tile_id)) deallocate(domain2d_spec_out%tile_id) !< Check if allocated
1866  allocate(domain2d_spec_out%tile_id(starting:ending))
1867  domain2d_spec_out%tile_id = domain2d_spec_in%tile_id
1868  endif
1869 
1870  if (associated(domain2d_spec_in%x)) then
1871  starting = lbound(domain2d_spec_in%x,1)
1872  ending = ubound(domain2d_spec_in%x,1)
1873 
1874  if (associated(domain2d_spec_out%x)) deallocate(domain2d_spec_out%x) !< Check if allocated
1875  allocate(domain2d_spec_out%x(starting:ending))
1876  do i = starting, ending
1877  call mpp_copy_domain1d_spec(domain2d_spec_in%x(i), domain2d_spec_out%x(i))
1878  enddo
1879  endif
1880 
1881  if (associated(domain2d_spec_in%y)) then
1882  starting = lbound(domain2d_spec_in%y,1)
1883  ending = ubound(domain2d_spec_in%y,1)
1884 
1885  if (associated(domain2d_spec_out%y)) deallocate(domain2d_spec_out%y) !< Check if allocated
1886  allocate(domain2d_spec_out%y(starting:ending))
1887  do i = starting, ending
1888  call mpp_copy_domain1d_spec(domain2d_spec_in%y(i), domain2d_spec_out%y(i))
1889  enddo
1890  endif
1891 
1892  end subroutine mpp_copy_domain2d_spec
1893 
1894  !> @brief Copies input 1d domain spec to the output 1d domain spec
1895  subroutine mpp_copy_domain1d_spec(domain1D_spec_in, domain1D_spec_out)
1896  type(domain1d_spec), intent(in) :: domain1D_spec_in !< Input
1897  type(domain1d_spec), intent(out) :: domain1D_spec_out !< Output
1898 
1899  domain1d_spec_out%pos = domain1d_spec_in%pos
1900 
1901  call mpp_copy_domain_axis_spec(domain1d_spec_in%compute, domain1d_spec_out%compute)
1902  call mpp_copy_domain_axis_spec(domain1d_spec_in%global, domain1d_spec_out%global)
1903  end subroutine mpp_copy_domain1d_spec
1904 
1905  !> @brief Copies input domain_axis_spec to the output domain_axis_spec
1906  subroutine mpp_copy_domain_axis_spec(domain_axis_spec_in, domain_axis_spec_out)
1907  type(domain_axis_spec), intent(in) :: domain_axis_spec_in !< Input
1908  type(domain_axis_spec), intent(out) :: domain_axis_spec_out !< Output
1909 
1910  domain_axis_spec_out%begin = domain_axis_spec_in%begin
1911  domain_axis_spec_out%end = domain_axis_spec_in%end
1912  domain_axis_spec_out%size = domain_axis_spec_in%size
1913  domain_axis_spec_out%max_size = domain_axis_spec_in%max_size
1914  domain_axis_spec_out%is_global = domain_axis_spec_in%is_global
1915  end subroutine mpp_copy_domain_axis_spec
1916 
1917  !######################################################################
1918  subroutine set_group_update(group, domain)
1919  type(mpp_group_update_type), intent(inout) :: group
1920  type(domain2d), intent(inout) :: domain
1921  integer :: nscalar, nvector, nlist
1922  integer :: nsend, nrecv, nsend_old, nrecv_old
1923  integer :: nsend_s, nsend_x, nsend_y
1924  integer :: nrecv_s, nrecv_x, nrecv_y
1925  integer :: update_buffer_pos, tot_recv_size, tot_send_size
1926  integer :: msgsize_s, msgsize_x, msgsize_y, msgsize
1927  logical :: recv_s(8), send_s(8)
1928  logical :: recv_x(8), send_x(8), recv_y(8), send_y(8)
1929  integer :: ntot, n, l, m, ksize
1930  integer :: i_s, i_x, i_y, rank_s, rank_x, rank_y, rank
1931  integer :: ind_s(3*MAXOVERLAP)
1932  integer :: ind_x(3*MAXOVERLAP)
1933  integer :: ind_y(3*MAXOVERLAP)
1934  integer :: pelist(3*MAXOVERLAP)
1935  integer :: send_size(3*MAXOVERLAP)
1936  integer :: position_x, position_y, npack, nunpack, dir
1937  integer :: pack_buffer_pos, unpack_buffer_pos
1938  integer :: omp_get_num_threads, nthreads
1939  character(len=8) :: text
1940  type(overlap_type), pointer :: overPtr => null()
1941  type(overlapspec), pointer :: update_s => null()
1942  type(overlapspec), pointer :: update_x => null()
1943  type(overlapspec), pointer :: update_y => null()
1944 
1945  nscalar = group%nscalar
1946  nvector = group%nvector
1947 
1948  !--- get the overlap data type
1949  select case(group%gridtype)
1950  case (agrid)
1951  position_x = center
1952  position_y = center
1953  case (bgrid_ne, bgrid_sw)
1954  position_x = corner
1955  position_y = corner
1956  case (cgrid_ne, cgrid_sw)
1957  position_x = east
1958  position_y = north
1959  case (dgrid_ne, dgrid_sw)
1960  position_x = north
1961  position_y = east
1962  case default
1963  call mpp_error(fatal, "set_group_update: invalid value of gridtype")
1964  end select
1965  if(nscalar>0) then
1966  update_s => search_update_overlap(domain, group%whalo_s, group%ehalo_s, &
1967  group%shalo_s, group%nhalo_s, group%position)
1968  endif
1969  if(nvector>0) then
1970  update_x => search_update_overlap(domain, group%whalo_v, group%ehalo_v, &
1971  group%shalo_v, group%nhalo_v, position_x)
1972  update_y => search_update_overlap(domain, group%whalo_v, group%ehalo_v, &
1973  group%shalo_v, group%nhalo_v, position_y)
1974  endif
1975 
1976  if(nscalar > 0) then
1977  recv_s = group%recv_s
1978  send_s = recv_s
1979  endif
1980  if(nvector > 0) then
1981  recv_x = group%recv_x
1982  send_x = recv_x
1983  recv_y = group%recv_y
1984  send_y = recv_y
1985  end if
1986  nlist = size(domain%list(:))
1987  group%initialized = .true.
1988  nsend_s = 0; nsend_x = 0; nsend_y = 0
1989  nrecv_s = 0; nrecv_x = 0; nrecv_y = 0
1990 
1991  if(nscalar > 0) then
1992  !--- This check could not be done because of memory domain
1993 ! if( group%isize_s .NE. (group%ie_s-group%is_s+1) .OR. group%jsize_s .NE. (group%je_s-group%js_s+1)) &
1994 ! call mpp_error(FATAL, "set_group_update: mismatch of size of the field and domain memory domain")
1995  nsend_s = update_s%nsend
1996  nrecv_s = update_s%nrecv
1997  endif
1998 
1999  !--- ksize_s must equal ksize_v
2000  if(nvector > 0 .AND. nscalar > 0) then
2001  if(group%ksize_s .NE. group%ksize_v) then
2002  call mpp_error(fatal, "set_group_update: ksize_s and ksize_v are not equal")
2003  endif
2004  ksize = group%ksize_s
2005  else if (nscalar > 0) then
2006  ksize = group%ksize_s
2007  else if (nvector > 0) then
2008  ksize = group%ksize_v
2009  else
2010  call mpp_error(fatal, "set_group_update: nscalar and nvector are all 0")
2011  endif
2012 
2013  nthreads = 1
2014 !$OMP PARALLEL
2015 !$ nthreads = omp_get_num_threads()
2016 !$OMP END PARALLEL
2017  if( nthreads > nthread_control_loop ) then
2018  group%k_loop_inside = .false.
2019  else
2020  group%k_loop_inside = .true.
2021  endif
2022 
2023  if(nvector > 0) then
2024  !--- This check could not be done because of memory domain
2025 ! if( group%isize_x .NE. (group%ie_x-group%is_x+1) .OR. group%jsize_x .NE. (group%je_x-group%js_x+1)) &
2026 ! call mpp_error(FATAL, "set_group_update: mismatch of size of the fieldx and domain memory domain")
2027 ! if( group%isize_y .NE. (group%ie_y-group%is_y+1) .OR. group%jsize_y .NE. (group%je_y-group%js_y+1)) &
2028 ! call mpp_error(FATAL, "set_group_update: mismatch of size of the fieldy and domain memory domain")
2029  nsend_x = update_x%nsend
2030  nrecv_x = update_x%nrecv
2031  nsend_y = update_y%nsend
2032  nrecv_y = update_y%nrecv
2033  endif
2034 
2035  !figure out message size for each processor.
2036  ntot = nrecv_s + nrecv_x + nrecv_y
2037  if(ntot > 3*maxoverlap) call mpp_error(fatal, "set_group_update: ntot is greater than 3*MAXOVERLAP")
2038  n = 1
2039  i_s = 1
2040  i_x = 1
2041  i_y = 1
2042  ind_s = -1
2043  ind_x = -1
2044  ind_y = -1
2045  nrecv = 0
2046  do while(n<=ntot)
2047  if( i_s <= nrecv_s ) then
2048  rank_s = update_s%recv(i_s)%pe-domain%pe
2049  if(rank_s .LE. 0) rank_s = rank_s + nlist
2050  else
2051  rank_s = -1
2052  endif
2053  if( i_x <= nrecv_x ) then
2054  rank_x = update_x%recv(i_x)%pe-domain%pe
2055  if(rank_x .LE. 0) rank_x = rank_x + nlist
2056  else
2057  rank_x = -1
2058  endif
2059  if( i_y <= nrecv_y ) then
2060  rank_y = update_y%recv(i_y)%pe-domain%pe
2061  if(rank_y .LE. 0) rank_y = rank_y + nlist
2062  else
2063  rank_y = -1
2064  endif
2065  nrecv = nrecv + 1
2066  rank = maxval((/rank_s, rank_x, rank_y/))
2067  if(rank == rank_s) then
2068  n = n + 1
2069  ind_s(nrecv) = i_s
2070  pelist(nrecv) = update_s%recv(i_s)%pe
2071  i_s = i_s + 1
2072  endif
2073  if(rank == rank_x) then
2074  n = n + 1
2075  ind_x(nrecv) = i_x
2076  pelist(nrecv) = update_x%recv(i_x)%pe
2077  i_x = i_x + 1
2078  endif
2079  if(rank == rank_y) then
2080  n = n + 1
2081  ind_y(nrecv) = i_y
2082  pelist(nrecv) = update_y%recv(i_y)%pe
2083  i_y = i_y + 1
2084  endif
2085  enddo
2086 
2087  nrecv_old = nrecv
2088  nrecv = 0
2089  update_buffer_pos = 0
2090  tot_recv_size = 0
2091 
2092  !--- setup for recv
2093  do l = 1, nrecv_old
2094  msgsize_s = 0
2095  msgsize_x = 0
2096  msgsize_y = 0
2097  m = ind_s(l)
2098  if(m>0) msgsize_s = get_mesgsize(update_s%recv(m), recv_s)*ksize*nscalar
2099  m = ind_x(l)
2100  if(m>0) msgsize_x = get_mesgsize(update_x%recv(m), recv_x)*ksize*nvector
2101  m = ind_y(l)
2102  if(m>0) msgsize_y = get_mesgsize(update_y%recv(m), recv_y)*ksize*nvector
2103  msgsize = msgsize_s + msgsize_x + msgsize_y
2104  if( msgsize.GT.0 )then
2105  tot_recv_size = tot_recv_size + msgsize
2106  nrecv = nrecv + 1
2107  if(nrecv > maxoverlap) then
2108  call mpp_error(fatal, "set_group_update: nrecv is greater than MAXOVERLAP, increase MAXOVERLAP")
2109  endif
2110  group%from_pe(nrecv) = pelist(l)
2111  group%recv_size(nrecv) = msgsize
2112  group%buffer_pos_recv(nrecv) = update_buffer_pos
2113  update_buffer_pos = update_buffer_pos + msgsize
2114  end if
2115  end do
2116  group%nrecv = nrecv
2117 
2118  !--- setup for unpack
2119  nunpack = 0
2120  unpack_buffer_pos = 0
2121  do l = 1, nrecv_old
2122  m = ind_s(l)
2123  if(m>0) then
2124  overptr => update_s%recv(m)
2125  do n = 1, overptr%count
2126  dir = overptr%dir(n)
2127  if(recv_s(dir)) then
2128  nunpack = nunpack + 1
2129  if(nunpack > maxoverlap) call mpp_error(fatal, &
2130  "set_group_update: nunpack is greater than MAXOVERLAP, increase MAXOVERLAP 1")
2131  group%unpack_type(nunpack) = field_s
2132  group%unpack_buffer_pos(nunpack) = unpack_buffer_pos
2133  group%unpack_rotation(nunpack) = overptr%rotation(n)
2134  group%unpack_is(nunpack) = overptr%is(n)
2135  group%unpack_ie(nunpack) = overptr%ie(n)
2136  group%unpack_js(nunpack) = overptr%js(n)
2137  group%unpack_je(nunpack) = overptr%je(n)
2138  group%unpack_size(nunpack) = overptr%msgsize(n)*nscalar
2139  unpack_buffer_pos = unpack_buffer_pos + group%unpack_size(nunpack)*ksize
2140  end if
2141  end do
2142  end if
2143 
2144  m = ind_x(l)
2145  if(m>0) then
2146  overptr => update_x%recv(m)
2147  do n = 1, overptr%count
2148  dir = overptr%dir(n)
2149  if(recv_x(dir)) then
2150  nunpack = nunpack + 1
2151  if(nunpack > maxoverlap) call mpp_error(fatal, &
2152  "set_group_update: nunpack is greater than MAXOVERLAP, increase MAXOVERLAP 2")
2153  group%unpack_type(nunpack) = field_x
2154  group%unpack_buffer_pos(nunpack) = unpack_buffer_pos
2155  group%unpack_rotation(nunpack) = overptr%rotation(n)
2156  group%unpack_is(nunpack) = overptr%is(n)
2157  group%unpack_ie(nunpack) = overptr%ie(n)
2158  group%unpack_js(nunpack) = overptr%js(n)
2159  group%unpack_je(nunpack) = overptr%je(n)
2160  group%unpack_size(nunpack) = overptr%msgsize(n)*nvector
2161  unpack_buffer_pos = unpack_buffer_pos + group%unpack_size(nunpack)*ksize
2162  end if
2163  end do
2164  end if
2165 
2166  m = ind_y(l)
2167  if(m>0) then
2168  overptr => update_y%recv(m)
2169  do n = 1, overptr%count
2170  dir = overptr%dir(n)
2171  if(recv_y(dir)) then
2172  nunpack = nunpack + 1
2173  if(nunpack > maxoverlap) call mpp_error(fatal, &
2174  "set_group_update: nunpack is greater than MAXOVERLAP, increase MAXOVERLAP 3")
2175  group%unpack_type(nunpack) = field_y
2176  group%unpack_buffer_pos(nunpack) = unpack_buffer_pos
2177  group%unpack_rotation(nunpack) = overptr%rotation(n)
2178  group%unpack_is(nunpack) = overptr%is(n)
2179  group%unpack_ie(nunpack) = overptr%ie(n)
2180  group%unpack_js(nunpack) = overptr%js(n)
2181  group%unpack_je(nunpack) = overptr%je(n)
2182  group%unpack_size(nunpack) = overptr%msgsize(n)*nvector
2183  unpack_buffer_pos = unpack_buffer_pos + group%unpack_size(nunpack)*ksize
2184  end if
2185  end do
2186  end if
2187  end do
2188  group%nunpack = nunpack
2189 
2190  if(update_buffer_pos .NE. unpack_buffer_pos ) call mpp_error(fatal, &
2191  .NE."set_group_update: update_buffer_pos unpack_buffer_pos")
2192 
2193  !figure out message size for each processor.
2194  ntot = nsend_s + nsend_x + nsend_y
2195  n = 1
2196  i_s = 1
2197  i_x = 1
2198  i_y = 1
2199  ind_s = -1
2200  ind_x = -1
2201  ind_y = -1
2202  nsend = 0
2203  do while(n<=ntot)
2204  if( i_s <= nsend_s ) then
2205  rank_s = update_s%send(i_s)%pe-domain%pe
2206  if(rank_s .LT. 0) rank_s = rank_s + nlist
2207  else
2208  rank_s = nlist+1
2209  endif
2210  if( i_x <= nsend_x ) then
2211  rank_x = update_x%send(i_x)%pe-domain%pe
2212  if(rank_x .LT. 0) rank_x = rank_x + nlist
2213  else
2214  rank_x = nlist+1
2215  endif
2216  if( i_y <= nsend_y ) then
2217  rank_y = update_y%send(i_y)%pe-domain%pe
2218  if(rank_y .LT. 0) rank_y = rank_y + nlist
2219  else
2220  rank_y = nlist+1
2221  endif
2222  nsend = nsend + 1
2223  rank = minval((/rank_s, rank_x, rank_y/))
2224  if(rank == rank_s) then
2225  n = n + 1
2226  ind_s(nsend) = i_s
2227  pelist(nsend) = update_s%send(i_s)%pe
2228  i_s = i_s + 1
2229  endif
2230  if(rank == rank_x) then
2231  n = n + 1
2232  ind_x(nsend) = i_x
2233  pelist(nsend) = update_x%send(i_x)%pe
2234  i_x = i_x + 1
2235  endif
2236  if(rank == rank_y) then
2237  n = n + 1
2238  ind_y(nsend) = i_y
2239  pelist(nsend) = update_y%send(i_y)%pe
2240  i_y = i_y + 1
2241  endif
2242  enddo
2243 
2244  nsend_old = nsend
2245  nsend = 0
2246  tot_send_size = 0
2247  do l = 1, nsend_old
2248  msgsize_s = 0
2249  msgsize_x = 0
2250  msgsize_y = 0
2251  m = ind_s(l)
2252  if(m>0) msgsize_s = get_mesgsize(update_s%send(m), send_s)*ksize*nscalar
2253  m = ind_x(l)
2254  if(m>0) msgsize_x = get_mesgsize(update_x%send(m), send_x)*ksize*nvector
2255  m = ind_y(l)
2256  if(m>0) msgsize_y = get_mesgsize(update_y%send(m), send_y)*ksize*nvector
2257  msgsize = msgsize_s + msgsize_x + msgsize_y
2258  if( msgsize.GT.0 )then
2259  tot_send_size = tot_send_size + msgsize
2260  nsend = nsend + 1
2261  if(nsend > maxoverlap) then
2262  call mpp_error(fatal, "set_group_update: nsend is greater than MAXOVERLAP, increase MAXOVERLAP")
2263  endif
2264  send_size(nsend) = msgsize
2265  group%to_pe(nsend) = pelist(l)
2266  group%buffer_pos_send(nsend) = update_buffer_pos
2267  group%send_size(nsend) = msgsize
2268  update_buffer_pos = update_buffer_pos + msgsize
2269  end if
2270  end do
2271  group%nsend = nsend
2272 
2273  !--- setup for pack
2274  npack = 0
2275  pack_buffer_pos = unpack_buffer_pos
2276  do l = 1, nsend_old
2277  m = ind_s(l)
2278  if(m>0) then
2279  overptr => update_s%send(m)
2280  do n = 1, overptr%count
2281  dir = overptr%dir(n)
2282  if(send_s(dir)) then
2283  npack = npack + 1
2284  if(npack > maxoverlap) call mpp_error(fatal, &
2285  "set_group_update: npack is greater than MAXOVERLAP, increase MAXOVERLAP 1")
2286  group%pack_type(npack) = field_s
2287  group%pack_buffer_pos(npack) = pack_buffer_pos
2288  group%pack_rotation(npack) = overptr%rotation(n)
2289  group%pack_is(npack) = overptr%is(n)
2290  group%pack_ie(npack) = overptr%ie(n)
2291  group%pack_js(npack) = overptr%js(n)
2292  group%pack_je(npack) = overptr%je(n)
2293  group%pack_size(npack) = overptr%msgsize(n)*nscalar
2294  pack_buffer_pos = pack_buffer_pos + group%pack_size(npack)*ksize
2295  end if
2296  end do
2297  end if
2298 
2299  m = ind_x(l)
2300  if(m>0) then
2301  overptr => update_x%send(m)
2302  do n = 1, overptr%count
2303  dir = overptr%dir(n)
2304  !--- nonsym_edge update is not for rotation of 90 or -90 degree ( cubic sphere grid )
2305  if( group%nonsym_edge .and. (overptr%rotation(n)==ninety .or. &
2306  overptr%rotation(n)==minus_ninety) ) then
2307  call mpp_error(fatal, 'set_group_update: flags=NONSYMEDGEUPDATE is not compatible '// &
2308  'with 90 or -90 degree rotation (normally cubic sphere grid' )
2309  endif
2310  if(send_x(dir)) then
2311  npack = npack + 1
2312  if(npack > maxoverlap) call mpp_error(fatal, &
2313  "set_group_update: npack is greater than MAXOVERLAP, increase MAXOVERLAP 2")
2314  group%pack_type(npack) = field_x
2315  group%pack_buffer_pos(npack) = pack_buffer_pos
2316  group%pack_rotation(npack) = overptr%rotation(n)
2317  group%pack_is(npack) = overptr%is(n)
2318  group%pack_ie(npack) = overptr%ie(n)
2319  group%pack_js(npack) = overptr%js(n)
2320  group%pack_je(npack) = overptr%je(n)
2321  group%pack_size(npack) = overptr%msgsize(n)*nvector
2322  pack_buffer_pos = pack_buffer_pos + group%pack_size(npack)*ksize
2323  end if
2324  end do
2325  end if
2326 
2327  m = ind_y(l)
2328  if(m>0) then
2329  overptr => update_y%send(m)
2330  do n = 1, overptr%count
2331  dir = overptr%dir(n)
2332  if( group%nonsym_edge .and. (overptr%rotation(n)==ninety .or. &
2333  overptr%rotation(n)==minus_ninety) ) then
2334  call mpp_error(fatal, 'set_group_update: flags=NONSYMEDGEUPDATE is not compatible '// &
2335  'with 90 or -90 degree rotation (normally cubic sphere grid' )
2336  endif
2337  if(send_y(dir)) then
2338  npack = npack + 1
2339  if(npack > maxoverlap) call mpp_error(fatal, &
2340  "set_group_update: npack is greater than MAXOVERLAP, increase MAXOVERLAP 3")
2341  group%pack_type(npack) = field_y
2342  group%pack_buffer_pos(npack) = pack_buffer_pos
2343  group%pack_rotation(npack) = overptr%rotation(n)
2344  group%pack_is(npack) = overptr%is(n)
2345  group%pack_ie(npack) = overptr%ie(n)
2346  group%pack_js(npack) = overptr%js(n)
2347  group%pack_je(npack) = overptr%je(n)
2348  group%pack_size(npack) = overptr%msgsize(n)*nvector
2349  pack_buffer_pos = pack_buffer_pos + group%pack_size(npack)*ksize
2350  end if
2351  end do
2352  end if
2353  end do
2354  group%npack = npack
2355  if(update_buffer_pos .NE. pack_buffer_pos ) call mpp_error(fatal, &
2356  .NE."set_group_update: update_buffer_pos pack_buffer_pos")
2357 
2358  !--- make sure the buffer is large enough
2359  mpp_domains_stack_hwm = max( mpp_domains_stack_hwm, tot_recv_size+tot_send_size )
2360 
2361  if( mpp_domains_stack_hwm.GT.mpp_domains_stack_size )then
2362  write( text,'(i8)' )mpp_domains_stack_hwm
2363  call mpp_error( fatal, 'set_group_update: mpp_domains_stack overflow, '// &
2364  'call mpp_domains_set_stack_size('//trim(text)//') from all PEs.' )
2365  end if
2366 
2367  group%tot_msgsize = tot_recv_size+tot_send_size
2368 
2369 end subroutine set_group_update
2370 
2371 
2372 !######################################################################
2373  subroutine mpp_clear_group_update(group)
2374  type(mpp_group_update_type), intent(inout) :: group
2375 
2376  group%nscalar = 0
2377  group%nvector = 0
2378  group%nsend = 0
2379  group%nrecv = 0
2380  group%npack = 0
2381  group%nunpack = 0
2382  group%initialized = .false.
2383 
2384  end subroutine mpp_clear_group_update
2385 
2386 !#####################################################################
2388  type(mpp_group_update_type), intent(in) :: group
2389  logical :: mpp_group_update_initialized
2390 
2391  mpp_group_update_initialized = group%initialized
2392 
2393  end function mpp_group_update_initialized
2394 
2395 !#####################################################################
2396  function mpp_group_update_is_set(group)
2397  type(mpp_group_update_type), intent(in) :: group
2398  logical :: mpp_group_update_is_set
2399 
2400  mpp_group_update_is_set = (group%nscalar > 0 .OR. group%nvector > 0)
2401 
2402  end function mpp_group_update_is_set
2403 !> @}
subroutine mpp_get_overlap(domain, action, p, is, ie, js, je, dir, rot, position)
Set user stack size.
subroutine mpp_get_neighbor_pe_2d(domain, direction, pe)
Return PE North/South/East/West of this PE-domain. direction must be NORTH, SOUTH,...
subroutine mpp_get_global_domains1d(domain, begin, end, size)
Set user stack size.
type(overlapspec) function, pointer search_bound_overlap(domain, position)
This routine finds the bound at certain position.
subroutine mpp_set_super_grid_indices(grid)
Modifies the indices in the domain_axis_spec type to those of the supergrid.
integer function mpp_get_domain_npes(domain)
Set user stack size.
logical function mpp_domain2d_eq(a, b)
Set user stack size.
integer function mpp_get_tile_npes(domain)
Returns number of processors used on current tile.
integer function get_rank_send(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
Set user stack size.
subroutine set_group_update(group, domain)
Set user stack size.
subroutine mpp_get_global_domain2d(domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, tile_count, position)
Set user stack size.
subroutine mpp_get_pelist1d(domain, pelist, pos)
Set user stack size.
subroutine mpp_set_global_domain1d(domain, begin, end, size)
Set user stack size.
subroutine mpp_get_layout2d(domain, layout)
Set user stack size.
subroutine mpp_get_tile_pelist(domain, pelist)
Get the processors list used on current tile.
subroutine mpp_get_domain_components(domain, x, y, tile_count)
Retrieve 1D components of 2D decomposition.
subroutine mpp_get_domain_extents1d(domain, xextent, yextent)
Set user stack size.
logical function mpp_group_update_is_set(group)
Set user stack size.
integer function get_vector_send(domain, update_x, update_y, ind_x, ind_y, start_pos, pelist)
Set user stack size.
subroutine mpp_copy_domain_axis_spec(domain_axis_spec_in, domain_axis_spec_out)
Copies input domain_axis_spec to the output domain_axis_spec.
integer function mpp_get_domain_tile_commid(domain)
Set user stack size.
subroutine mpp_set_global_domain2d(domain, xbegin, xend, ybegin, yend, xsize, ysize, tile_count)
Set user stack size.
logical function mpp_domain1d_eq(a, b)
Set user stack size.
integer function get_rank_recv(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
Set user stack size.
integer function, dimension(size(domain%tile_id(:))) mpp_get_tile_id(domain)
Returns the tile_id on current pe.
logical function mpp_domain_is_symmetry(domain)
Set user stack size.
subroutine mpp_create_super_grid_domain(domain)
Modifies the indices of the input domain to create the supergrid domain.
integer function mpp_get_num_overlap(domain, action, p, position)
Set user stack size.
type(overlapspec) function, pointer search_check_overlap(domain, position)
this routine finds the check at certain position
integer function mpp_get_current_ntile(domain)
Returns number of tile on current pe.
subroutine mpp_set_domain_symmetry(domain, symmetry)
Set user stack size.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
subroutine mpp_get_domain_extents2d(domain, xextent, yextent)
This will return xextent and yextent for each tile.
subroutine mpp_get_global_domain1d(domain, begin, end, size, max_size)
Set user stack size.
logical function domain_update_is_needed(domain, whalo, ehalo, shalo, nhalo)
Set user stack size.
integer function, dimension(2) mpp_get_io_domain_layout(domain)
Set user stack size.
subroutine mpp_get_neighbor_pe_1d(domain, direction, pe)
Return PE to the righ/left of this PE-domain.
subroutine mpp_get_tile_compute_domains(domain, xbegin, xend, ybegin, yend, position)
Set user stack size.
subroutine mpp_get_compute_domain2d(domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, x_is_global, y_is_global, tile_count, position)
Set user stack size.
subroutine mpp_get_compute_domain1d(domain, begin, end, size, max_size, is_global)
Set user stack size.
subroutine mpp_copy_domain1d_spec(domain1D_spec_in, domain1D_spec_out)
Copies input 1d domain spec to the output 1d domain spec.
logical function mpp_domain1d_ne(a, b)
Set user stack size.
subroutine mpp_get_layout1d(domain, layout)
Set user stack size.
subroutine mpp_copy_domain2d_spec(domain2D_spec_in, domain2d_spec_out)
Copies input 2d domain spec to the output 2d domain spec.
subroutine mpp_get_update_pelist(domain, action, pelist, position)
Set user stack size.
subroutine mpp_set_data_domain2d(domain, xbegin, xend, ybegin, yend, xsize, ysize, x_is_global, y_is_global, tile_count)
Set user stack size.
type(overlapspec) function, pointer search_update_overlap(domain, whalo, ehalo, shalo, nhalo, position)
this routine found the domain has the same halo size with the input whalo, ehalo,
subroutine mpp_get_compute_domains1d(domain, begin, end, size)
Set user stack size.
integer function mpp_get_domain_tile_root_pe(domain)
Set user stack size.
subroutine mpp_get_pelist2d(domain, pelist, pos)
Set user stack size.
subroutine mpp_get_memory_domain2d(domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, x_is_global, y_is_global, position)
Set user stack size.
integer function get_mesgsize(overlap, do_dir)
Set user stack size.
subroutine mpp_set_compute_domain1d(domain, begin, end, size, is_global)
Set user stack size.
integer function mpp_get_domain_pe(domain)
Set user stack size.
integer function get_vector_recv(domain, update_x, update_y, ind_x, ind_y, start_pos, pelist)
Set user stack size.
subroutine mpp_get_tile_list(domain, tiles)
Return the tile_id on current pelist. one-tile-per-pe is assumed.
logical function mpp_group_update_initialized(group)
Set user stack size.
subroutine mpp_set_compute_domain2d(domain, xbegin, xend, ybegin, yend, xsize, ysize, x_is_global, y_is_global, tile_count)
Set user stack size.
type(domain2d) function, pointer mpp_get_io_domain(domain)
Set user stack size.
subroutine mpp_get_compute_domains2d(domain, xbegin, xend, xsize, ybegin, yend, ysize, position)
Set user stack size.
subroutine compute_overlaps(domain, position, update, check, ishift, jshift, x_cyclic_offset, y_cyclic_offset, whalo, ehalo, shalo, nhalo)
Computes remote domain overlaps.
subroutine mpp_copy_domain2d(domain_in, domain_out)
Copies input 2d domain to the output 2d domain.
subroutine mpp_get_memory_domain1d(domain, begin, end, size, max_size, is_global)
Set user stack size.
subroutine mpp_get_domain_shift(domain, ishift, jshift, position)
Returns the shift value in x and y-direction according to domain position..
logical function mpp_domain2d_ne(a, b)
Set user stack size.
logical function mpp_domain_is_initialized(domain)
Set user stack size.
subroutine nullify_domain2d_list(domain)
Set user stack size.
logical function mpp_domain_is_tile_root_pe(domain)
Returns if current pe is the root pe of the tile, if number of tiles on current pe is greater than 1,...
character(len=name_length) function mpp_get_domain_name(domain)
Set user stack size.
subroutine mpp_get_global_domains2d(domain, xbegin, xend, xsize, ybegin, yend, ysize, position)
Set user stack size.
integer function mpp_get_domain_commid(domain)
Set user stack size.
subroutine mpp_get_data_domain2d(domain, xbegin, xend, ybegin, yend, xsize, xmax_size, ysize, ymax_size, x_is_global, y_is_global, tile_count, position)
Set user stack size.
subroutine mpp_get_domain_pelist(domain, pelist)
Set user stack size.
subroutine mpp_domains_set_stack_size(n)
Set user stack size.
subroutine mpp_clear_group_update(group)
Set user stack size.
subroutine mpp_get_update_size(domain, nsend, nrecv, position)
Set user stack size.
integer function get_rank_unpack(domain, overlap_x, overlap_y, rank_x, rank_y, ind_x, ind_y)
Set user stack size.
integer function mpp_get_domain_root_pe(domain)
Set user stack size.
subroutine mpp_set_data_domain1d(domain, begin, end, size, is_global)
Set user stack size.
recursive subroutine mpp_copy_domain1d(domain_in, domain_out)
Copies input 1d domain to the output 1d domain.
subroutine set_overlaps(domain, overlap_in, overlap_out, whalo_out, ehalo_out, shalo_out, nhalo_out)
this routine sets up the overlapping for mpp_update_domains for arbitrary halo update....
subroutine mpp_get_data_domain1d(domain, begin, end, size, max_size, is_global)
Set user stack size.
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406