Skip to content
Snippets Groups Projects
Commit 23845225 authored by Thomas Jahns's avatar Thomas Jahns :cartwheel:
Browse files

Add distributed curvilinear grids to the test suite.

parent 20b1c11d
No related branches found
No related tags found
No related merge requests found
......@@ -63,11 +63,6 @@ modelRegionCompute(double region[], int nlev, int nlat, int nlon,
mscale, mrscale);
}
#ifdef USE_MPI
static void
findPartition2D(int npart[2], int num_parts);
#endif
void
modelRun(struct model_config setup, MPI_Comm comm)
{
......@@ -436,65 +431,6 @@ modelRun(struct model_config setup, MPI_Comm comm)
Free(filename);
}
#ifdef USE_MPI
static void
findPartition2D(int npart[2], int num_parts)
{
const uint64_t rscale = 256;
uint32_t *factors = NULL;
xassert(num_parts > 0);
int numFactors
= PPM_prime_factorization_32((uint32_t)num_parts, &factors);
/* try to distribute prime factors on dimensions to get
* approx. 2 times as many parts in x dim than y dim */
const uint64_t optimumRatio = rscale * 2;
npart[0] = num_parts, npart[1] = 1;
uint_fast32_t npart_attempt[2];
uint64_t bestRatio = (uint64_t)num_parts * rscale,
bestDiff = (uint64_t)llabs((long long)(bestRatio - optimumRatio));
/* test all assignments of factors to dimensions, starting with
* only one assigned to x dim (omitting 0 because that would
* always give npart[1] > npart[0] */
for (int assign2X = 1; assign2X <= numFactors; ++assign2X)
{
uint_fast32_t pattern = (UINT32_C(1) << assign2X) - 1,
lastPattern = pattern << (numFactors - assign2X);
do {
npart_attempt[0] = 1;
npart_attempt[1] = 1;
/* loop over all factors */
for (uint_fast32_t i = 0; i < (uint_fast32_t)numFactors; ++i)
{
uint_fast32_t dim_idx = (pattern >> i) & 1;
npart_attempt[dim_idx] *= factors[i];
}
uint64_t ratio = ((uint64_t)npart_attempt[0] * rscale)
/ (uint64_t)npart_attempt[1];
uint64_t diff = (uint64_t)llabs((long long)(ratio - optimumRatio));
if (diff < bestDiff)
{
npart[0] = (int)npart_attempt[0];
npart[1] = (int)npart_attempt[1];
bestDiff = diff;
bestRatio = ratio;
}
{
uint_fast32_t t;
#if HAVE_DECL___BUILTIN_CTZ
t = pattern | (pattern - 1);
pattern = (t + 1)
| (((~t & -~t) - 1) >> (__builtin_ctz((unsigned)pattern) + 1));
#else
t = (pattern | (pattern - 1)) + 1;
pattern = t | ((((t & -t) / (pattern & -pattern)) >> 1) - 1);
#endif
}
} while (pattern <= lastPattern);
}
Free(factors);
}
#endif
/*
* Local Variables:
* c-file-style: "Java"
......
......@@ -105,7 +105,7 @@ struct boolOptionParse
static struct boolOptionParse
parseBooleanLongOption(size_t optionStrSize, const char *optionStr,
const char *argStr)
const char *argStr)
{
bool invert;
struct boolOptionParse result;
......@@ -187,7 +187,7 @@ parse_long_option(struct model_config *restrict setup,
useDistGridOptionStr, str)).matched)
{
#if defined USE_MPI && defined HAVE_PPM_DIST_ARRAY_H
setup->flags = (setup->flags & ~~PIO_WRITE_CONFIG_USE_DIST_GRID_FLAG)
setup->flags = (setup->flags & ~PIO_WRITE_CONFIG_USE_DIST_GRID_FLAG)
| (bop.value << PIO_WRITE_CONFIG_USE_DIST_GRID_BIT);
#else
invalidOptionDie("CDI-PIO option -q%s feature unavailable %s\n",
......
#! @SHELL@
set -e
LOG="${LOG-pio_write.log}"
LOG=${LOG-pio_write.log}
if [ @ENABLE_GRIB@ != yes ]; then
# skip tests for unsupported formats
exit 77
......@@ -10,11 +10,13 @@ suffix="${suffix-grb}"
if [ "@USE_MPI@" = yes ]; then
variations="-qcache-redists -qcreate-curvilinear-grid"
else
unset variations
variations="-qcreate-curvilinear-grid"
fi
for variation in "" $variations ; do
test_variation()
{
exec 5>&1 6>&2 >"$LOG" 2>&1
echo "creating data" >&2
variation=$1
../libtool --mode=execute \
@MPI_LAUNCH@ \
-n ${mpi_task_num} ${tool_wrap} ./pio_write ${pio_write_args} ${variation}
......@@ -24,9 +26,18 @@ for variation in "" $variations ; do
echo "checking example_1.$suffix" >&2
../libtool --mode=execute \
${tool_wrap} ./cksum_read example_1.${suffix} example_1.cksum
exec 2>&6 1>&5 5>&- 6>&-
\rm "$LOG"
}
exec 5>&1 6>&2 >"$LOG" 2>&1
for variation in "" $variations ; do
test_variation "$variation"
done
if [ x"@USE_MPI@" = xyes -a x"@HAVE_PPM_DIST_ARRAY@" = xyes ]; then
test_variation "-qcreate-curvilinear-grid -quse-dist-grid"
fi
exec 2>&6 1>&5 5>&- 6>&-
\rm "$LOG"
#
# Local Variables:
# mode: sh
......
......@@ -2,35 +2,196 @@
# include "config.h"
#endif
#include <stdio.h>
#include <stdlib.h>
#ifdef USE_MPI
#include <mpi.h>
#include <yaxt.h>
#include "error.h"
#include "pio_util.h"
#endif
#include "dmemory.h"
#ifdef HAVE_PPM_CORE
#include <ppm/ppm_uniform_partition.h>
#include <core/ppm_combinatorics.h>
#endif
#include "cdi.h"
#ifdef USE_MPI
#include "cdipio.h"
#endif
#include "pio_write_setup_grid.h"
#include "cdi_uuid.h"
#include "simple_model_helper.h"
#if USE_MPI
void
findPartition2D(int npart[2], int num_parts)
{
const uint64_t rscale = 256;
uint32_t *factors = NULL;
xassert(num_parts > 0);
int numFactors
= PPM_prime_factorization_32((uint32_t)num_parts, &factors);
/* try to distribute prime factors on dimensions to get
* approx. 2 times as many parts in x dim than y dim */
const uint64_t optimumRatio = rscale * 2;
npart[0] = num_parts, npart[1] = 1;
uint_fast32_t npart_attempt[2];
uint64_t bestRatio = (uint64_t)num_parts * rscale,
bestDiff = (uint64_t)llabs((long long)(bestRatio - optimumRatio));
/* test all assignments of factors to dimensions, starting with
* only one assigned to x dim (omitting 0 because that would
* always give npart[1] > npart[0] */
for (int assign2X = 1; assign2X <= numFactors; ++assign2X)
{
uint_fast32_t pattern = (UINT32_C(1) << assign2X) - 1,
lastPattern = pattern << (numFactors - assign2X);
do {
npart_attempt[0] = 1;
npart_attempt[1] = 1;
/* loop over all factors */
for (uint_fast32_t i = 0; i < (uint_fast32_t)numFactors; ++i)
{
uint_fast32_t dim_idx = (pattern >> i) & 1;
npart_attempt[dim_idx] *= factors[i];
}
uint64_t ratio = ((uint64_t)npart_attempt[0] * rscale)
/ (uint64_t)npart_attempt[1];
uint64_t diff = (uint64_t)llabs((long long)(ratio - optimumRatio));
if (diff < bestDiff)
{
npart[0] = (int)npart_attempt[0];
npart[1] = (int)npart_attempt[1];
bestDiff = diff;
bestRatio = ratio;
}
{
uint_fast32_t t;
#if HAVE_DECL___BUILTIN_CTZ
t = pattern | (pattern - 1);
pattern = (t + 1)
| (((~t & -~t) - 1) >> (__builtin_ctz((unsigned)pattern) + 1));
#else
t = (pattern | (pattern - 1)) + 1;
pattern = t | ((((t & -t) / (pattern & -pattern)) >> 1) - 1);
#endif
}
} while (pattern <= lastPattern);
}
Free(factors);
}
#endif
int
setupGrid(struct model_config *setup, MPI_Comm comm)
{
int nlon = setup->nlon, nlat = setup->nlat;
int gridID;
if (setup->flags & PIO_WRITE_CONFIG_CREATE_CURVILINEAR_GRID_FLAG)
gridID = createLocalCurvilinearGrid(nlon, nlat);
else
gridID = createGlobalLatLonGrid(nlon, nlat);
if (setup->flags & PIO_WRITE_CONFIG_CREATE_UUID_FLAG)
{
int rank = 0;
int rank = 0;
#if USE_MPI
xmpi(MPI_Comm_rank(comm, &rank));
#else
(void)comm;
#endif
if (setup->flags & PIO_WRITE_CONFIG_USE_DIST_GRID_FLAG)
{
#if USE_MPI
xmpi(MPI_Comm_rank(comm, &rank));
if (setup->flags & PIO_WRITE_CONFIG_CREATE_CURVILINEAR_GRID_FLAG)
{
size_t gridsize = (size_t)nlon * (size_t)nlat;
Xt_int grid_global_shape[2] = { (Xt_int)nlat, (Xt_int)nlon },
grid_local_chunk_start[2];
int nprocs;
xmpi(MPI_Comm_size(comm, &nprocs));
int nprocs2D[2], grid_local_chunk_shape[2];
findPartition2D(nprocs2D, nprocs);
{
struct PPM_extent set_interval_y = { 0, nlat },
part_interval_y
= PPM_uniform_partition(set_interval_y, nprocs2D[1],
rank / nprocs2D[0]);
grid_local_chunk_shape[0] = part_interval_y.size;
grid_local_chunk_start[0] = part_interval_y.first;
}
{
struct PPM_extent set_interval_x = { 0, nlon },
part_interval_x
= PPM_uniform_partition(set_interval_x, nprocs2D[0],
rank % nprocs2D[0]);
grid_local_chunk_shape[1] = part_interval_x.size;
grid_local_chunk_start[1] = part_interval_x.first;
}
{
Xt_idxlist partDesc2D
= xt_idxsection_new((Xt_int)0, 2, grid_global_shape,
grid_local_chunk_shape,
grid_local_chunk_start);
const int gridDecoChunks[2][2] = { [0] = {
[0] = (int)grid_local_chunk_start[0],
[1] = (int)grid_local_chunk_shape[0],
}, [1] = {
[0] = (int)grid_local_chunk_start[1],
[1] = (int)grid_local_chunk_shape[1],
} };
gridID = cdiPioDistGridCreate(GRID_CURVILINEAR, (int)gridsize,
nlon, nlat, 0, gridDecoChunks,
partDesc2D, NULL, NULL);
xt_idxlist_delete(partDesc2D);
}
size_t grid_chunk_size = (size_t)grid_local_chunk_shape[0]
* (size_t)grid_local_chunk_shape[1];
{
double *grid_chunk = (double *)
Malloc(2 * grid_chunk_size * sizeof (*grid_chunk));
/* anti-clockwise coordinates around Amazonia */
static const struct cart_coord region[4]
#ifdef __cplusplus
= { { DEG2RAD(-25.0), DEG2RAD(-85.0) },
{ DEG2RAD(-18.0), DEG2RAD(-44.0) },
{ DEG2RAD(-7.0), DEG2RAD(-50.0) },
{ DEG2RAD(10.0), DEG2RAD(-80.0) } };
#else
= { { .lon = DEG2RAD(-85.0), .lat = DEG2RAD(-25.0) },
{ .lon = DEG2RAD(-44.0), .lat = DEG2RAD(-18.0) },
{ .lon = DEG2RAD(-50.0), .lat = DEG2RAD(7.0) },
{ .lon = DEG2RAD(-80.0), .lat = DEG2RAD(10.0) } };
#endif
const size_t xyRange[2][2]
= { { (size_t)grid_local_chunk_start[1],
(size_t)grid_local_chunk_start[0] },
{ (size_t)grid_local_chunk_shape[1],
(size_t)grid_local_chunk_shape[0] } };
computeCurvilinearChunk(grid_chunk, region,
(size_t)nlon, (size_t)nlat, xyRange);
gridDefXvals(gridID, grid_chunk + grid_chunk_size);
gridDefYvals(gridID, grid_chunk);
Free(grid_chunk);
}
}
else
{
fputs("Error: distributed lat-lon grids are unsupported"
" at this moment.\n", stderr);
abort();
}
#else
(void)comm;
fputs("Error: distributed grids are unsupported for"
" non-MPI configurations.\n", stderr);
abort();
#endif
}
else
{
if (setup->flags & PIO_WRITE_CONFIG_CREATE_CURVILINEAR_GRID_FLAG)
gridID = createLocalCurvilinearGrid(nlon, nlat);
else
gridID = createGlobalLatLonGrid(nlon, nlat);
}
if (setup->flags & PIO_WRITE_CONFIG_CREATE_UUID_FLAG)
{
unsigned char uuid[CDI_UUID_SIZE];
if (rank == 0)
cdiCreateUUID(uuid);
......@@ -42,6 +203,8 @@ setupGrid(struct model_config *setup, MPI_Comm comm)
return gridID;
}
/*
* Local Variables:
* c-file-style: "Java"
......
......@@ -11,6 +11,11 @@
#include "pio_write.h"
#ifdef USE_MPI
void
findPartition2D(int npart[2], int num_parts);
#endif
int
setupGrid(struct model_config *setup, MPI_Comm comm);
......
......@@ -9,8 +9,14 @@
#include <stdlib.h>
#include <string.h>
#ifdef HAVE_PPM_CORE
#include <ppm/ppm_uniform_partition.h>
#include <core/ppm_combinatorics.h>
#endif
#include "cdi.h"
#include "dmemory.h"
#include "error.h"
#include "simple_model_helper.h"
......@@ -140,16 +146,6 @@ createGlobalLatLonGrid(int nlon, int nlat)
return gridID;
}
struct cart_coord {
double lat, lon;
};
static void
computeCurvilinearChunk(double *coords_,
const struct cart_coord a[4],
size_t sizex, size_t sizey,
const size_t xyRange[2][2]);
int
createLocalCurvilinearGrid(int sizex, int sizey)
{
......@@ -213,7 +209,7 @@ intermediateCoord(struct cart_coord p1, struct cart_coord p2, double f)
return ic;
}
static void
void
computeCurvilinearChunk(double *coords_,
const struct cart_coord a[4],
size_t sizex, size_t sizey,
......@@ -221,10 +217,10 @@ computeCurvilinearChunk(double *coords_,
{
size_t startx = xyRange[0][0],
starty = xyRange[0][1],
endx = xyRange[1][0],
endy = xyRange[1][1],
chunkSizeX = endx - startx,
chunkSizeY = endy - starty;
chunkSizeX = xyRange[1][0],
chunkSizeY = xyRange[1][1],
endx = startx + chunkSizeX,
endy = starty + chunkSizeY;
assert(startx <= endx && endx <= sizex && starty <= endy && endy <= sizey);
#ifdef __cplusplus
auto coords = (double (*)[chunkSizeY][chunkSizeX])coords_;
......@@ -247,9 +243,7 @@ computeCurvilinearChunk(double *coords_,
}
}
#if defined USE_MPI && ! defined HAVE_PPM_CORE
static int32_t
uniform_partition_start(struct PPM_extent set_interval, int nparts,
int part_idx);
......@@ -302,7 +296,6 @@ PPM_prime_factorization_32(uint32_t n, uint32_t **factors)
*factors = pfactors;
return numFactors;
}
#endif
/*
......
......@@ -57,7 +57,28 @@ createGlobalLatLonGrid(int nlon, int nlat);
int
createLocalCurvilinearGrid(int sizex, int sizey);
#if defined (USE_MPI) && ! defined(HAVE_PPM_CORE)
struct cart_coord {
double lat, lon;
};
/*
* coords_: points to array of size xyRange[0][0] * xyRange[0][1] * 2, where
* the first half of xyRange[0][0] * xyRange[0][1] elements will be
* initialized to longitude values and the second half to the latitude values
* a: cartesian coordinates of rectangular part of globe in
* anti-clockwise order
* sizex, sizey: number of grid points in x and y directions respectively
* xyRange: subset of the ranges [0,sizex) and [0,sizey) to compute,
* given as start in xyRange[0] and count values in xyRange[1] for x
* at index 0 and y at index 1
*/
void
computeCurvilinearChunk(double *coords_,
const struct cart_coord a[4],
size_t sizex, size_t sizey,
const size_t xyRange[2][2]);
#if defined USE_MPI && ! defined HAVE_PPM_CORE
struct PPM_extent
{
int32_t first, size;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment