Contact: B. Wyman
Reviewers:
Change history: WebCVS Log for bgrid_horiz.f90
Initializes horizontal grid constants needed by the B-grid dynamical core
and determines the domain decomposition for running on distributed
memory machines. This routine calls the FMS mpp_domains package.
A data structure (derived-type variable, horiz_grid_type) is returned
that contains the grid constants and mpp domain data structures.
As it's name implies, the model uses the staggered Arakawa B-grid.
The model's prognostic variables are carried on two separate sub-grids.
The momentum components are located together on the momentum or
velocity grid; while pressure, temperature, and tracers are located
together on the mass or temperature grid. These two sub-grid will
use the module's public parameters VGRID and TGRID, respectively.
More detail about the the B-grid is found in the notes section.
mpp_domains_mod
constants_mod
fms_mod
use bgrid_horiz_mod [ ,only: horiz_grid_type,
bgrid_type,
horiz_grid_init,
get_horiz_grid_size,
get_horiz_grid_bound,
update_np, update_sp,
TGRID, VGRID ]
bgrid_type, horiz_grid_type
Data structues (i.e., derived-type variables) containing
horizontal grid constants needed by the B-grid dynamical core.
horiz_grid_init
Function that initializes/returns a variable of
type (horiz_grid_type).
get_horiz_grid_size
Subroutine that returns the number of grid boxes along the x and
y axes of the temperature grid or velocity grid. There is an
option for returning either the global or local compute grid sizes.
Usually used in conjunction with routine get_horiz_grid_bound.
get_horiz_grid_bound
Subroutine that returns the longitude and latitude boundaries (edges)
of temperature or velocity grid boxes. There is an option for
returning either the global or local compute grid edges.
Routine get_horiz_grid_size is usually called first to get
the axis sizes.
update_np, update_sp
Logical functions that determine if the northern (or southern)
halo rows for the local processor lie within the global halo
rows (ie.e, need a polar boundary update).
TGRID, VGRID
Integer parameters to be used as the "grid" argument with interfaces
get_horiz_grid_size, get_horiz_grid_bound, update_np, and update_sp.
is, ie starting, ending x-axis index in the compute domain [integer,scalar]
js, je starting, ending y-axis index in the compute domain [integer,scalar]
isd, ied starting, ending x-axis index in the data domain [integer,scalar]
jsd, jed starting, ending y-axis index in the data domain [integer,scalar]
isg, ieg starting, ending x-axis index in the global domain [integer,scalar]
jsg, jeg starting, ending y-axis index in the global domain [integer,scalar]
dx local grid spacing for x-axis (in meters) [real,dimension(ilb:iub,jlb:jub)]
rdx reciprocal of dx (1/m) [real,dimension(jlb:jub)]
dy local grid spacing for y-axis (in meters) [real,scalar]
rdy reciprocal of dy (1/m) [real,scalar]
area area of a local grid box (in m2) [real,dimension(jlb:jub)]
rarea reciprocal of area (1/m2) [real,dimension(jlb:jub)]
areasum bit-reproducible global sum of area (in m2) [real,scalar]
tlm longitude at the center of a local grid box (in radians)
[real,dimension(ilb:iub)]
tph latitude at the center of a local grid box (in radians)
[real,dimension(jlb:jub)]
aph, alm actual latitude, longitude at the center of a local grid box (in radians)
[real,dimension(ilb:iub,jlb:jub)]
blatg latitude boundaries of grid boxes along the global y-axis (in radians)
[real,dimension(jsg:jeg+1)]
blong longitude boundaries grid boxes along the global x-axis (in radians)
[real,dimension(isg:ieg+1)]
Domain domain2D derived-type variable with halosize halo size = 1
[type(domain2d)]
Domain_nohalo domain2D derived-type variable with halo size = 0, used for outputting diagnostic fields
[type(domain2d)]
Tmp grid constants for the temperature/tracer/mass grid [type(bgrid_type)]
Vel grid constants for the u/v wind component grid [type(bgrid_type)]
Interp interpolation weights used by module bgrid_change_grid [type(bgrid_interp_type)]
sinphv sine of Vel%aph [real,dimension(ilb:iub,jlb:jub)]
tanphv tangent of Vel%aph [real,dimension(ilb:iub,jlb:jub)]
nlon number of longitude (x-axis) grid points in the global grid (no halo points)
[integer,scalar]
nlat number of latitude (y-axis) grid points in the global grid (no halo points)
[integer,scalar]
isize, jsize size of arrays on the local processor's grid (including halo points)
Note: isize=iub-ilb+1, jsize=jub-jlb+1
[integer,scalar]
ilb, iub lower, upper bounds of local x-axis indexing [integer,scalar]
jlb, jub lower, upper bounds of local y-axis indexing [integer,scalar]
ihalo, jhalo number of halo points along the east and west boundaries (ihalo) or
south and north boundaries (jhalo) [integer,scalar]
dlm grid spacing for x-axis (in radians) [real,scalar]
dph grid spacing for y-axis (in radians) [real,scalar]
dlmd grid spacing for x-axis (in degrees of longitude) [real,scalar]
dphd grid spacing for y-axis (in degrees of latitude) [real,scalar]
sb southern boundary of temperature grid (in radians) [real,scalar]
nb northern boundary of temperature grid (in radians) [real,scalar]
wb western boundary of temperature grid (in radians) [real,scalar]
eb eastern boundary of temperature grid (in radians) [real,scalar]
decompx indicates if the x-axis has been decomposed across processors [logical,scalar]
decompy indicates if the y-axis has been decomposed across processors [logical,scalar]
channel indicates if the channel model feature has been implemented (internal option)
[logical,scalar]
double_periodic indicates if the double periodic (f-plane) channel feature has been select (internal option)
[logical,scalar]
NOTES: The local indexing variables: is, ie, js, je, isd, ied, jsd, jed, ilb, iub, jlb, jub, are consistent with the global index values (isg, ieg, jsg, jeg). All local arrays (on either the temperature or velocity grid) have the same horizontal size [dimension(ilb:iub,jlb:jub)]. Due to the grid staggering, the northernmost velocity domains do not use the last j-row (jub). The domain boundaries: sb, nb, wb, eb, are only used when channel = true. Otherwise the global domain boundaries: sb=-pi/2, nb=pi/2, wb=0, eb=2pi are used.
call horiz_grid_init ( Hgrid, nlon, nlat [, decomp] )
INPUT
nlon, nlat The number of global longitude, latitude grid points (respectively)
for the mass/temperature grid.
[integer,scalar]
INPUT/OUTPUT
Hgrid Data structure that contains all necessary horizontal grid information
and domain decomposition variables needed by the model.
[type(horiz_grid_type)]
OPTIONAL INPUT
decomp The domain decomposition for the x and y axis, where
decomp(1) = x-axis decomposition, decomp(2) = y-axis
decomposition. If decomp(1)*decomp(2) does not equal
the number of processors the model will fail.
If either or both decomp(1), decomp(2) is zero,
the rules below apply.
[integer,dimension(2), default: decomp=0,0]
NOTES
1) If decomp(1)=0 and decomp(2)>0, then decomp(1)=decomp(2)/(# PEs), or vice versa.
The final decomposition must always satisfy: decomp(1)*decomp(2)=number of processors.
2) If decomp(1)=decomp(2)=0, then MPP_DEFINE_LAYOUT will be called.
It is likely that two-dimensional decomposition will be implemented.
call get_horiz_grid_size ( Hgrid, grid, nlon, nlat [, global] )
INPUT
Hgrid Data structure returned by a previous call to horiz_grid_init.
[type(horiz_grid_type)]
grid Specifies which grid (temperature or velocity) the returned grid size
will be for. Use the publicly accessible parameters: TGRID or VGRID.
[integer,scalar]
OUTPUT
nlon, nlat The number of grid points along the x- and y-axis, respectively.
The returned values will be for either the global grid size or
the local processor's grid size depending on the value of optional
argument "global".
[integer,scalar]
OPTIONAL INPUT
global Flag that determines if the returned values are for the global
grid (global=.TRUE.) or the local compute grid (global=.FALSE.).
[logical,scalar, default: FALSE]
call get_horiz_grid_bound ( Hgrid, grid, blon, blat [, global] )
INPUT
Hgrid Data structure returned by a previous call to horiz_grid_init.
[type(horiz_grid_type)]
grid Specifies which grid (temperature or velocity) the returned grid boundaries
will be for. Use the publicly accessible parameters: TGRID or VGRID.
[integer,scalar]
OUTPUT
blon Longitude edges of grid boxes for either the global grid or the
local processor's grid depending on the value of optional argument "global".
[real,dimension(nlon+1)]
blat Latitude edges of grid boxes for either the global grid or the
local processor's grid depending on the value of optional argument "global".
[real,dimension(nlat+1)]
OPTIONAL INPUT
global Flag that determines if the returned values are for the global
grid (global=.TRUE.) or the local compute grid (global=.FALSE.).
[logical,scalar, default: FALSE]
answer = update_np ( Hgrid, grid )
INPUT
Hgrid The data structure returned by a previous call to horiz_grid_init.
[type(horiz_grid_type)]
grid Specifies which grid (temperature or velocity) the returned value
will be for. Use the publicly accessible parameters: TGRID or VGRID.
[integer,scalar]
RETURNED VALUE
Returns TRUE when the halo rows along the north boundary lie within
global halo rows, otherwise FALSE is returned. Note that for single
processor runs the returned value will always be TRUE.
[logical]
answer = update_sp ( Hgrid, grid )
INPUT
Hgrid The data structure returned by a previous call to horiz_grid_init.
[type(horiz_grid_type)]
grid Specifies which grid (temperature or velocity) the returned value
will be for. Use the publicly accessible parameters: TGRID or VGRID.
[integer,scalar]
RETURNED VALUE
Returns TRUE when the halo rows along the south boundary lie within
global halo rows, otherwise FALSE is returned. Note that for single
processor runs the returned value will always be TRUE.
[logical]
FATAL errors in horiz_grid_init
number of processor requested not compatible with grid
The domain decomposition either requested or computed was
not compatible with the resolution of the grid.
If you requested a decomposition other than the default,
check to make sure it is correct.
Also, the number of processor you are running on may not be
compatible with the resolution of the model.
Otherwise there may be a code error.
FATAL errors in get_horiz_grid_bound
invalid argument dimension for blon
The size of the blon output argument must equal the
number of longitude grid boxes plus one.
invalid argument dimension for blat
The size of the blat output argument must equal the
number of latitude grid boxes plus one.
FATAL errors in get_horiz_grid_size, get_horiz_grid_bound,
update_sp, update_np
invalid grid
The input argument "grid" has an incorrect value.
Make sure you are using one of public parameters, TGRID or VGRID.
Arakawa, A. and V. R. Lamb, 1977: Computational design of the basic dynamical processes of the UCLA general circulation model. Methods in Computational Physics, Vol. 17. Academic Press, pp. 173-265.
* The grid transformation option (tph0d, tlm0d) has not been
extensively tested. If want to use this option check with
the developer.
* The channel model mode has not been extensively tested.
THE B-GRID
----------
The B-grid can be visualized as two overlapping grids, one that
contains momentum (u and v) and the other mass fields (surface
pressure, temperature, and tracers). These separate grids are
rectangular in shape, defined by longitudes along the x-axis and
latitude along the y-axis. The grids are diagonally shifted from
each other, such that, the center of a momentum grid box is located
at the corners where four mass grid boxes intersect.
Horizontal indexing increases from west to east, and from south to north.
Indexing is set up so that a velocity point with the same i,j is
located to the north and east of the mass (temperature) point.
T(i-1,j+1) T(i,j+1) T(i+1,j+1)
v(i-1,j ) v(i,j )
T(i-1,j ) T(i,j ) T(i+1,j )
v(i-1,j-1) v(i,j-1)
T(i-1,j-1) T(i,j-1) T(i+1,j-1)
For computational purposes, extra rows and columns of grid boxes
are carried on the perimeter of the horizontal grid. These extra
points are called halo points. The number of halo points along
the west/east boundaries and south/north boundaries are specified
as separate arguments when initializing the horizontal grid. The
default is for one halo row along all boundaries. Because of the
grid staggering, the global momentum grid does not use it's
northernmost row.
The interior grid (excluding halo points) is called the compute domain.
The global size of the longitude axis (first dimension) can have an
even or odd number of points. The latitude axis (2nd dimension) should
always be an even number of points. When there is an even number of
global latitude points, mass points will straddle the equator,
while velocity points lie along the equator.
UNSUPPORTED INTERNAL OPTIONS
----------------------------
tph0d,tlm0d Latitude/longitude for shifting the position of the poles.
Set both to zero for no transformation (the default).
This option has not been recently used and is currently
not supported.
[real, default: tph0d=0., tlm0d=0.]
wbd western edge in degrees of the first longitude (i=1) of temperature grid boxes
[real, default: wbd = 0.0]
ebd eastern edge in degrees of the last longitude (i=nlon) of temperature grid boxes
[real, default: ebd = 360.0]
sbd southern edge in degrees of the first latitude row (j=1) of temperature grid boxes
[real, default: sbd = -90.0]
nbd northern edge in degrees of the last latitude row (j=nlat) of temperature grid boxes
[real, default: nbd = 90.0]
do_channel if TRUE, then use sbd,nbd as southern and northern boundaries of a channel,
and wbd,ebd as the periodic western and eastern boundaries.
[logical, default: do_channel=FALSE]
do_fplane_approx f-plane channel without spherical geometry used in conjunction
with option do_channel = .TRUE.
[logical, default: do_fplane_approx=FALSE]
fphd coriolis latitude in degrees used for f-plane channel
[real, default: fphd=45.]
do_double_periodic domain is cyclic in X and Y, this option must be used with
do_fplane_approx = .TRUE. Tmp and Vel grids will have the same
size global compute domain.
[logical, default: do_double_periodic=FALSE]
None.