mpp_domains_mod is a set of simple calls for domain decomposition and domain updates on rectilinear grids. It requires the module , upon which it is built.
type, public :: domain_axis_spec private integer :: begin, end, size, max_size logical :: is_global end type domain_axis_spec type, public :: domain1D private type(domain_axis_spec) :: compute, data, global, active logical :: mustputb, mustgetb, mustputf, mustgetf, folded type(domain1D), pointer, dimension(:) :: list integer :: pe !PE to which this domain is assigned integer :: pos end type domain1D domaintypes of higher rank can be constructed from type domain1D typically we only need 1 and 2D, but could need higher (e.g 3D LES) some elements are repeated below if they are needed once per domain type, public :: domain2D private type(domain1D) :: x type(domain1D) :: y type(domain2D), pointer, dimension(:) :: list integer :: pe !PE to which this domain is assigned integer :: pos end type domain2D type(domain1D), public :: NULL_DOMAIN1D type(domain2D), public :: NULL_DOMAIN2DThe domain2D type contains all the necessary information to define the global, compute and data domains of each task, as well as the PE associated with the task. The PEs from which remote data may be acquired to update the data domain are also contained in a linked list of neighbours.
mpp_mod
use mpp_domains_mod [, only: mpp_define_domains,
mpp_update_domains,
mpp_redistribute,
mpp_global_field,
mpp_global_max,
mpp_global_sum,
operator,
mpp_get_compute_domain,
mpp_get_compute_domains,
mpp_get_data_domain,
mpp_get_global_domain,
mpp_define_layout,
mpp_get_pelist,
mpp_get_layout,
mpp_domains_init,
mpp_domains_set_stack_size,
mpp_domains_exit,
mpp_get_domain_components ]
call mpp_define_domains ( global_indices, ndivs, domain, & pelist, flags, halo, extent, maskmap )
call mpp_define_domains ( global_indices, layout, domain, pelist, & xflags, yflags, xhalo, yhalo, & xextent, yextent, maskmap, name )
global_indices | Defines the global domain. [integer, dimension(2)] [integer, dimension(4)] |
ndivs | Is the number of domain divisions required. [integer] |
pelist | List of PEs to which the domains are to be assigned. [integer, dimension(0:)] [integer, dimension(0:)] |
flags | An optional flag to pass additional information
about the desired domain topology. Useful flags in a 1D decomposition
include GLOBAL_DATA_DOMAIN and CYCLIC_GLOBAL_DOMAIN. Flags are integers: multiple flags may
be added together. The flag values are public parameters available by
use association. [integer] |
halo | Width of the halo. [integer] |
extent | Normally mpp_define_domains attempts
an even division of the global domain across ndivs domains. The extent array can be used by the user to pass a
custom domain division. The extent array has ndivs elements and holds the compute domain widths, which should add up to
cover the global domain exactly. [integer, dimension(0:)] |
maskmap | Some divisions may be masked
(maskmap=.FALSE.) to exclude them from the computation (e.g
for ocean model domains that are all land). The maskmap array
is dimensioned ndivs and contains .TRUE. values for
any domain that must be included in the computation (default
all). The pelist array length should match the number of
domains included in the computation. [logical, dimension(0:)] [logical, dimension(:,:)] |
layout |
[integer, dimension(2)] |
xflags, yflags |
[integer] |
xhalo, yhalo |
[integer] |
xextent, yextent |
[integer, dimension(0:)] |
name |
[character(len=*)] |
domain | Holds the resulting domain decomposition. [type(domain1D)] [type(domain2D)] |
call mpp_define_domains( (/1,100/), 10, domain, & flags=GLOBAL_DATA_DOMAIN+CYCLIC_GLOBAL_DOMAIN, halo=2 )defines 10 compute domains spanning the range [1,100] of the global domain. The compute domains are non-overlapping blocks of 10. All the data domains are global, and with a halo of 2 span the range [-1:102]. And since the global domain has been declared to be cyclic, domain(9)%next => domain(0) and domain(0)%prev => domain(9). A field is allocated on the data domain, and computations proceed on the compute domain. A call to would fill in the values in the halo region:
call mpp_get_data_domain( domain, isd, ied ) !returns -1 and 102 call mpp_get_compute_domain( domain, is, ie ) !returns (1,10) on PE 0 ... allocate( a(isd:ied) ) do i = is,ie a(i) = <perform computations> end do call mpp_update_domains( a, domain )The call to mpp_update_domains fills in the regions outside the compute domain. Since the global domain is cyclic, the values at i=(-1,0) are the same as at i=(99,100); and i=(101,102) are the same as i=(1,2).
call mpp_define_domains( (/1,100,1,100/), (/2,2/), domain, xhalo=1 )will create the following domain layout:
|---------|-----------|-----------|-------------| |domain(1)|domain(2) |domain(3) |domain(4) | |--------------|---------|-----------|-----------|-------------| |Compute domain|1,50,1,50|51,100,1,50|1,50,51,100|51,100,51,100| |--------------|---------|-----------|-----------|-------------| |Data domain |0,51,1,50|50,101,1,50|0,51,51,100|50,101,51,100| |--------------|---------|-----------|-----------|-------------|Again, we allocate arrays on the data domain, perform computations on the compute domain, and call mpp_update_domains to update the halo region.
call mpp_define_domains( (/1,100,1,100/), layout=(/4,1/), domain, xhalo=1 )This will create the following domain layout:
|----------|-----------|-----------|------------| |domain(1) |domain(2) |domain(3) |domain(4) | |--------------|----------|-----------|-----------|------------| |Compute domain|1,100,1,25|1,100,26,50|1,100,51,75|1,100,76,100| |--------------|----------|-----------|-----------|------------| |Data domain |0,101,1,25|0,101,26,50|0,101,51,75|1,101,76,100| |--------------|----------|-----------|-----------|------------|
call mpp_update_domains ( field, domain, flags )
call mpp_update_domains ( fieldx, fieldy, domain, flags, gridtype )
call mpp_redistribute ( domain_in, field_in, domain_out, field_out )
field_in | field_in is dimensioned on the data domain of domain_in. |
field_out | field_out on the data domain of domain_out. |
call mpp_global_field ( domain, local, global, flags )
domain | |
local | local is dimensioned on either the compute domain or the data domain of domain. |
flags | flags can be given the value XONLY or YONLY, to specify a globalization on one axis only. |
global | global is dimensioned on the corresponding global domain. |
mpp_global_max ( domain, field, locus )
domain | |
field | field is dimensioned on either the compute domain or the data domain of domain. |
locus | locus, if present, can be used to retrieve the location of the maximum (as in the MAXLOC intrinsic of f90). |
call mpp_global_sum ( domain, field, flags )
type(domain1D) :: a, b type(domain2D) :: c, d ... if( a.NE.b )then ... end if if( c==d )then ... end ifDomains are considered equal if and only if the start and end indices of each of their component global, data and compute domains are equal.
call mpp_get_compute_domain
call mpp_get_compute_domains ( domain, xbegin, xend, xsize, & ybegin, yend, ysize )
domain |
xbegin,ybegin | |
xend,yend | |
xsize,ysize |
call mpp_get_data_domain
call mpp_get_global_domain
call mpp_define_layout ( global_indices, ndivs, layout )
global_indices |
[integer, dimension(4)] |
ndivs |
[integer] |
layout |
[integer, dimension(2)] |
domain |
[type(domain1D)] [type(domain2D)] |
pelist |
[integer, dimension(:)] [integer, dimension(:)] |
pos |
[integer] [integer] |
call mpp_get_layout ( domain, layout )
domain |
[type(domain1D)] [type(domain2D)] |
layout |
[integer] [integer, dimension(2)] |
call mpp_domains_init (flags)
flags |
[integer] |
call mpp_domains_set_stack_size (n)
n |
[integer] |
call mpp_domains_exit ()
call mpp_get_domain_components ( domain, x, y )
domain |
[type(domain2D)] |
x,y |
[type(domain1D)] |
use mpp_domains_modmpp_domains_mod uses , and therefore is subject to the compiling and linking requirements of that module.
The source consists of the main source file and also requires the following include files: GFDL users can check it out of the main CVS repository as part of the CVS module. The current public tag is . External users can download the latest package . Public access to the GFDL CVS repository will soon be made available.