A test for the parallel file interpolation method.
#include <stdlib.h>
#include <unistd.h>
#include "tests.h"
#include "test_common.h"
#include "dist_grid_utils.h"
#include "weight_file_common.h"
#include <mpi.h>
#include <yaxt.h>
#include <netcdf.h>
double const err_tol = 1e-14;
char const src_grid_name[] = "src_grid";
char const tgt_grid_name[] = "tgt_grid";
char const file_name[] = "test_interp_method_file_parallel_weights.nc";
char const file_name_2[] = "test_interp_method_file_parallel_weights_2.nc";
static void target_main(MPI_Comm target_comm);
static void source_main(MPI_Comm source_comm);
setenv("YAC_IO_MAX_NUM_RANKS_PER_NODE", "2", 1);
MPI_Init(NULL, NULL);
xt_initialize(MPI_COMM_WORLD);
int comm_rank, comm_size;
MPI_Comm_rank(MPI_COMM_WORLD, &comm_rank);
MPI_Comm_size(MPI_COMM_WORLD, &comm_size);
MPI_Barrier(MPI_COMM_WORLD);
if (comm_size != 3) {
PUT_ERR("ERROR: wrong number of processes");
xt_finalize();
MPI_Finalize();
return TEST_EXIT_CODE;
}
int tgt_flag = comm_rank < 1;
MPI_Comm split_comm;
MPI_Comm_split(MPI_COMM_WORLD, tgt_flag, 0, &split_comm);
if (tgt_flag) target_main(split_comm);
else source_main(split_comm);
MPI_Comm_free(&split_comm);
xt_finalize();
MPI_Finalize();
return TEST_EXIT_CODE;
}
static void source_main(MPI_Comm source_comm) {
int my_source_rank;
MPI_Comm_rank(source_comm, &my_source_rank);
{
if (my_source_rank == 0) {
int src_indices[] = {0,1,2,3,4,5,6,7,8,9,10,11};
int tgt_indices[] = {0,1,2,3,4,5,6,7,8,9,10,11};
double weights[] = {0,1,2,3,4,5,6,7,8,9,10,11};
size_t num_links = 12;
unsigned num_src_fields = 2;
int num_links_per_field[2] = {num_links, 0};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int with_halo = 1;
int global_corner_mask[3][4] = {
{1,1,1,1}, {1,1,0,0}, {0,0,0,0}};
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
int * src_corner_mask =
src_corner_mask[i] =
((
int*)(&(global_corner_mask[0][0])))[grid_data.
vertex_ids[i]];
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
{
++collection_idx) {
src_data[collection_idx] =
xmalloc(1 *
sizeof(**src_data));
src_data[collection_idx][0] =
src_data[collection_idx][0][i] =
(double)(collection_idx * 12);
}
++collection_idx) {
free(src_data[collection_idx][0]);
free(src_data[collection_idx]);
}
free(src_data);
}
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int tgt_indices[] = { 0, 0, 0, 0,
1, 1, 1, 1,
2, 2, 2, 2,
3, 3, 3, 3,
4, 4, 4, 4,
5, 5, 5, 5};
int src_indices[] = { 0, 1, 4, 5,
1, 2, 5, 6,
2, 3, 6, 7,
4, 5, 8, 9,
5, 6, 9,10,
6, 7,10,11};
double weights[] = {0.1,0.2,0.3,0.4,
0.5,0.6,0.7,0.8,
0.9,1.0,1.1,1.2,
1.3,1.4,1.5,1.6,
1.7,1.8,1.9,2.0,
2.1,2.2,2.3,2.4};
size_t num_links = 24;
unsigned num_src_fields = 2;
int num_links_per_field[2] = {num_links, 0};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
{
double ** src_fields =
xmalloc(1 *
sizeof(*src_fields));
double *** src_data = &src_fields;
double * src_field =
src_fields[0] = src_field;
free(src_fields[0]);
free(src_fields);
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int tgt_indices[] = { 0, 0, 0, 0,
1, 1, 1, 1,
2, 2, 2, 2,
3, 3, 3, 3,
4, 4, 4, 4,
5, 5, 5, 5};
int src_indices[] = { 0, 1, 4, 5,
1, 2, 5, 6,
2, 3, 6, 7,
4, 5, 8, 9,
5, 6, 9,10,
6, 7,10,11};
double weights[] = {0.1,0.2,0.3,0.4,
0.5,0.6,0.7,0.8,
0.9,1.0,1.1,1.2,
1.3,1.4,1.5,1.6,
1.7,1.8,1.9,2.0,
2.1,2.2,2.3,2.4};
size_t num_links = 24;
unsigned num_src_fields = 2;
int num_links_per_field[2] = {num_links, 0};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int global_src_vertex_mask[] = {1, 1, 1, 0,
1, 1, 1, 1,
0, 1, 1, 1};
int global_src_cell_mask[] = {1, 1, 0,
0, 1, 1};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
int * src_vertex_mask =
xmalloc(12 *
sizeof(*src_vertex_mask));
int * src_cell_mask =
xmalloc(6 *
sizeof(*src_cell_mask));;
src_vertex_mask[i] =
for (
size_t i = 0; i < grid_data.
num_cells; ++i)
src_cell_mask[i] =
global_src_cell_mask[grid_data.
cell_ids[i]];
{{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = 0},
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = 0}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
{
double ** src_fields =
xmalloc(1 *
sizeof(*src_fields));
double *** src_data = &src_fields;
double * src_field =
src_fields[0] = src_field;
free(src_fields[0]);
free(src_fields);
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int * tgt_indices = NULL;
int * src_indices = NULL;
double * weights = NULL;
size_t num_links = 0;
unsigned num_src_fields = 2;
int num_links_per_field[2] = {num_links, 0};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
{
double ** src_fields =
xmalloc(1 *
sizeof(*src_fields));
double *** src_data = &src_fields;
double * src_field =
src_fields[0] = src_field;
free(src_fields[0]);
free(src_fields);
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int * tgt_indices = NULL;
int * src_indices = NULL;
double * weights = NULL;
size_t num_links = 0;
unsigned num_src_fields = 2;
int num_links_per_field[2] = {num_links, 0};
int tgt_id_fixed[] = {0, 2, 4, 6};
size_t num_fixed_tgt = 4;
double fixed_values[] = {-1.0, -2.0};
int num_tgt_per_fixed_value[] = {1, 3};
size_t num_fixed_values = 2;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
{
double ** src_fields =
xmalloc(1 *
sizeof(*src_fields));
double *** src_data = &src_fields;
double * src_field =
src_fields[0] = src_field;
free(src_fields[0]);
free(src_fields);
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int * tgt_indices = NULL;
int * src_indices = NULL;
double * weights = NULL;
size_t num_links = 0;
unsigned num_src_fields = 1;
int num_links_per_field[1] = {num_links};
int tgt_id_fixed[] = {1,3,5,7,9,11,0,2,4,6,8,10};
size_t num_fixed_tgt = 12;
double fixed_values[] = {-1.0, -2.0};
int num_tgt_per_fixed_value[] = {6, 6};
size_t num_fixed_values = 2;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {3,2};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,2},{1,2}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
{
double ** src_fields =
xmalloc(1 *
sizeof(*src_fields));
double *** src_data = &src_fields;
double * src_field =
src_fields[0] = src_field;
free(src_fields[0]);
free(src_fields);
}
}
if (my_source_rank == 0) unlink(file_name);
}
{
if (my_source_rank == 0) {
int src_indices[3][36] = {{0,1,2,3,
0,1,3,4, 1,2,4,5, 3,4,6,7, 4,5,7,8,
0,1,3,5, 2,3,4,7, 5,6,8,10, 7,8,9,11},
{0,1,2,3,
0,1,3,4, 1,2,4,5, 3,4,6,7, 4,5,7,8,
0,1,3,5, 2,3,4,7, 5,6,8,10, 7,8,9,11},
{1,2,
4,
11}};
int tgt_indices[3][36] = {{0,1,2,3,
0,0,0,0, 1,1,1,1, 2,2,2,2, 3,3,3,3,
0,0,0,0, 1,1,1,1, 2,2,2,2, 3,3,3,3},
{0,1,2,3,
0,0,0,0, 1,1,1,1, 2,2,2,2, 3,3,3,3,
0,0,0,0, 1,1,1,1, 2,2,2,2, 3,3,3,3},
{1,2,
0,
3}};
double weights[3][36] = {{1,1,1,1,
0.25,0.25,0.25,0.25, 0.25,0.25,0.25,0.25,
0.25,0.25,0.25,0.25, 0.25,0.25,0.25,0.25,
0.25,0.25,0.25,0.25, 0.25,0.25,0.25,0.25,
0.25,0.25,0.25,0.25, 0.25,0.25,0.25,0.25},
{1,1,1,1,
1,1,1,1, 1,1,1,1,
1,1,1,1, 1,1,1,1,
1,1,1,1, 1,1,1,1,
1,1,1,1, 1,1,1,1},
{1,1,1,1}};
size_t num_links[3] = {1*4 + 4*4 + 4*4,
1*4 + 4*4 + 4*4,
1*2 + 1*1 + 1*1};
unsigned num_src_fields = 3;
int num_links_per_field[3][3] = {{1*4, 4*4, 4*4},
{1*4, 4*4, 4*4},
{2, 1, 1}};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
num_src_fields, num_links_per_field[
weight_type], tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0};
double coordinates_y[] = {0.0, 1.0, 2.0};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2][2] = {{0,0},{1,0}};
size_t local_count[2][2] = {{1,2},{1,2}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_EDGE, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
weights, file_name_2, src_grid_name, tgt_grid_name, 0, 0);
for (int from_file = 0; from_file <= 1; ++from_file) {
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
(from_file)?weights_from_file:weights,
{
++collection_idx) {
double ** src_fields =
xmalloc(3 *
sizeof(*src_fields));
src_data[collection_idx] = src_fields;
{
double * src_field =
for (
size_t i = 0; i < grid_data.
num_cells; ++i)
src_field[i] =
(
double)(grid_data.
cell_ids[i] + 10 * collection_idx) + 0.0;
src_fields[0] = src_field;
}
{
double * src_field =
src_field[i] =
(
double)(grid_data.
vertex_ids[i] + 10 * collection_idx) + 0.3;
src_fields[1] = src_field;
}
{
double * src_field =
for (
size_t i = 0; i < grid_data.
num_edges; ++i)
src_field[i] =
(
double)(grid_data.
edge_ids[i] + 10 * collection_idx) + 0.7;
src_fields[2] = src_field;
}
}
++collection_idx) {
for (size_t i = 0; i < 3; ++i) free(src_data[collection_idx][i]);
free(src_data[collection_idx]);
}
free(src_data);
}
}
}
}
if (my_source_rank == 0) {
unlink(file_name);
unlink(file_name_2);
}
}
}
{
if (my_source_rank == 0) {
int src_indices[16] = {0,1,4,8,13,14,
8,12,13,17,24,
8,17,25,26,30};
int tgt_indices[16] = {0,1,4,8,13,14,
2,5,6,9,15,
3,7,10,11,12};
double weights[16] = {1,1,1,1,1,1,
1,1,1,1,1,
1,1,1,1,1};
size_t num_links = 16;
unsigned num_src_fields = 3;
int num_links_per_field[3] = {6, 5, 5};
int * tgt_id_fixed = NULL;
size_t num_fixed_tgt = 0;
double * fixed_values = NULL;
int * num_tgt_per_fixed_value = NULL;
size_t num_fixed_values = 0;
write_weight_file(
file_name, src_indices, tgt_indices, weights, num_links,
src_locations, num_src_fields, num_links_per_field, tgt_id_fixed,
num_fixed_tgt, fixed_values, num_tgt_per_fixed_value,
num_fixed_values, tgt_location, src_grid_name, tgt_grid_name);
}
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0, 4.0};
double coordinates_y[] = {0.0, 1.0, 2.0, 3.0, 4.0};
size_t const num_global_cells[2] = {4,4};
size_t local_start[2][2] = {{0,0},{2,0}};
size_t local_count[2][2] = {{2,4},{2,4}};
int with_halo = 1;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start[my_source_rank], local_count[my_source_rank], with_halo);
{{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_EDGE, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
weights, file_name_2, src_grid_name, tgt_grid_name, 0, 0);
for (int from_file = 0; from_file <= 1; ++from_file) {
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
(from_file)?weights_from_file:weights,
{
++collection_idx) {
double ** src_fields =
xmalloc(3 *
sizeof(*src_fields));
src_data[collection_idx] = src_fields;
{
double * src_field =
for (
size_t i = 0; i < grid_data.
num_cells; ++i)
src_field[i] =
(
double)(grid_data.
cell_ids[i] + 10 * collection_idx) + 0.0;
src_fields[0] = src_field;
}
{
double * src_field =
src_field[i] =
(
double)(grid_data.
vertex_ids[i] + 10 * collection_idx) + 0.3;
src_fields[1] = src_field;
}
{
double * src_field =
for (
size_t i = 0; i < grid_data.
num_edges; ++i)
src_field[i] =
(
double)(grid_data.
edge_ids[i] + 10 * collection_idx) + 0.7;
src_fields[2] = src_field;
}
}
++collection_idx) {
for (size_t i = 0; i < 3; ++i) free(src_data[collection_idx][i]);
free(src_data[collection_idx]);
}
free(src_data);
}
}
}
}
if (my_source_rank == 0) {
unlink(file_name);
unlink(file_name_2);
}
}
}
static void target_main(MPI_Comm target_comm) {
int my_target_rank;
MPI_Comm_rank(target_comm, &my_target_rank);
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = 0},
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
++collection_idx) {
tgt_data[collection_idx] =
tgt_data[collection_idx][i] = -1;
}
++collection_idx)
if (i >= 6) {
if (tgt_data[collection_idx][i] != -1.0)
PUT_ERR("wrong interpolation result");
} else if (fabs((double)(i * (i + 12 * collection_idx)) -
tgt_data[collection_idx][i]) > 1e-9) {
PUT_ERR("wrong interpolation result");
}
++collection_idx)
free(tgt_data[collection_idx]);
free(tgt_data);
}
}
}
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
double ** tgt_data =
xmalloc(1 *
sizeof(*tgt_data));
tgt_data[0] =
tgt_data[0][i] = -1;
double ref_tgt_field[9] =
{0.1*0+0.2*1+0.3*4+0.4*5,
0.5*1+0.6*2+0.7*5+0.8*6,
0.9*2+1.0*3+1.1*6+1.2*7,
1.3*4+1.4*5+1.5*8+1.6*9,
1.7*5+1.8*6+1.9*9+2.0*10,
2.1*6+2.2*7+2.3*10+2.4*11,
-1,-1,-1};
if (i >= 6) {
if (tgt_data[0][i] != -1.0)
PUT_ERR("wrong interpolation result");
} else if (fabs(ref_tgt_field[i] - tgt_data[0][i]) > 1e-9) {
PUT_ERR("wrong interpolation result");
}
free(tgt_data[0]);
free(tgt_data);
}
}
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = 0}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
double ** tgt_data =
xmalloc(1 *
sizeof(*tgt_data));
tgt_data[0] =
tgt_data[0][i] = -1;
double ref_tgt_field[9] =
{0.1*0+0.2*1+0.3*4+0.4*5,
0.5*1+0.6*2+0.7*5+0.8*6,
-1,
-1,
1.7*5+1.8*6+1.9*9+2.0*10,
2.1*6+2.2*7+2.3*10+2.4*11,
-1,-1,-1};
if (i >= 6) {
if (tgt_data[0][i] != -1.0)
PUT_ERR("wrong interpolation result");
} else if (fabs(ref_tgt_field[i] - tgt_data[0][i]) > 1e-9) {
PUT_ERR("wrong interpolation result");
}
free(tgt_data[0]);
free(tgt_data);
}
}
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
double ** tgt_data =
xmalloc(1 *
sizeof(*tgt_data));
tgt_data[0] =
tgt_data[0][i] = -1;
if (tgt_data[0][i] != -1.0)
PUT_ERR("wrong interpolation result");
free(tgt_data[0]);
free(tgt_data);
}
}
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CELL, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
double ** tgt_data =
xmalloc(1 *
sizeof(*tgt_data));
tgt_data[0] =
for (
size_t i = 0; i < grid_data.
num_vertices; ++i) tgt_data[0][i] = -3;
double ref_tgt_field[9] =
{-1,-3,-2, -3,-2,-3, -2,-3,-3};
if (ref_tgt_field[i] != tgt_data[0][i])
PUT_ERR("wrong interpolation result");
free(tgt_data[0]);
free(tgt_data);
}
}
{
double coordinates_x[] = {0.5, 1.5, 2.5};
double coordinates_y[] = {0.5, 1.5, 2.5};
size_t const num_global_cells[2] = {2,2};
size_t local_start[2] = {0,0};
size_t local_count[2] = {2,2};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
weights, reorder_type[i], 1,
double ** tgt_data =
xmalloc(1 *
sizeof(*tgt_data));
for (
size_t i = 0; i < grid_data.
num_edges; ++i) tgt_data[0][i] = -3;
double ref_tgt_field[12] = {-2,-1,-2,-1,-2,-1,-2,-1,-2,-1,-2,-1};
for (
size_t i = 0; i < grid_data.
num_edges; ++i)
if (ref_tgt_field[i] != tgt_data[0][i])
PUT_ERR("wrong interpolation result");
free(tgt_data[0]);
free(tgt_data);
}
}
{
double coordinates_x[] = {0.5, 1.5};
double coordinates_y[] = {0.5, 1.5};
size_t const num_global_cells[2] = {1,1};
size_t local_start[2] = {0,0};
size_t local_count[2] = {1,1};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_EDGE, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
weights, file_name_2, src_grid_name, tgt_grid_name, 0, 0);
for (int from_file = 0; from_file <= 1; ++from_file) {
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
(from_file)?weights_from_file:weights,
++collection_idx) {
tgt_data[collection_idx] =
tgt_data[collection_idx][i] = -3;
}
double ref_tgt_field[3][4] =
{{0.0 + 0.25 * (0.3+1.3+3.3+4.3) + 0.25 * (0.7+1.7+3.7+ 5.7),
1.0 + 0.25 * (1.3+2.3+4.3+5.3) + 0.25 * (2.7+3.7+4.7+ 7.7),
2.0 + 0.25 * (3.3+4.3+6.3+7.3) + 0.25 * (5.7+6.7+8.7+10.7),
3.0 + 0.25 * (4.3+5.3+7.3+8.3) + 0.25 * (7.7+8.7+9.7+11.7)},
{0.0 + (0.3+1.3+3.3+4.3) + (0.7+1.7+3.7+ 5.7),
1.0 + (1.3+2.3+4.3+5.3) + (2.7+3.7+4.7+ 7.7),
2.0 + (3.3+4.3+6.3+7.3) + (5.7+6.7+8.7+10.7),
3.0 + (4.3+5.3+7.3+8.3) + (7.7+8.7+9.7+11.7)},
{4.3, 1.0, 2.0, 11.7}};
double collection_factor[3] = {30.0, 90.0, 10.0};
collection_idx++)
(double)collection_idx) -
tgt_data[collection_idx][i]) > 1e-9)
PUT_ERR("wrong interpolation result");
++collection_idx) free(tgt_data[collection_idx]);
free(tgt_data);
}
}
}
}
}
{
double coordinates_x[] = {0.0, 1.0, 2.0, 3.0, 4.0};
double coordinates_y[] = {0.0, 1.0, 2.0, 3.0, 4.0};
size_t const num_global_cells[2] = {4,4};
size_t local_start[2] = {0,0};
size_t local_count[2] = {4,4};
int with_halo = 0;
for (size_t i = 0; i <= num_global_cells[0]; ++i)
for (size_t i = 0; i <= num_global_cells[1]; ++i)
yac_generate_basic_grid_data_reg2d(
coordinates_x, coordinates_y, num_global_cells,
local_start, local_count, with_halo);
{.location =
YAC_LOC_CORNER, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX},
{.location =
YAC_LOC_EDGE, .coordinates_idx = SIZE_MAX, .masks_idx = SIZE_MAX}};
size_t num_src_fields = sizeof(src_fields) / sizeof(src_fields[0]);
weights, file_name_2, src_grid_name, tgt_grid_name, 0, 0);
for (int from_file = 0; from_file <= 1; ++from_file) {
for (size_t i = 0; i < sizeof(reorder_type) / sizeof(reorder_type[0]);
++i) {
(from_file)?weights_from_file:weights,
++collection_idx) {
tgt_data[collection_idx] =
for (
size_t i = 0; i < grid_data.
num_cells; ++i)
tgt_data[collection_idx][i] = -3;
}
double ref_tgt_field[16] =
{ 0.0, 1.0, 8.3, 8.7,
4.0, 12.3, 13.3, 17.7,
8.0, 17.3, 25.7, 26.7,
30.7, 13.0, 14.0, 24.3};
collection_idx++)
for (
size_t i = 0; i < grid_data.
num_cells; ++i)
if (fabs((ref_tgt_field[i] + 10.0 * (double)collection_idx) -
tgt_data[collection_idx][i]) > 1e-9)
PUT_ERR("wrong interpolation result");
++collection_idx) free(tgt_data[collection_idx]);
free(tgt_data);
}
}
}
}
}
struct yac_basic_grid * yac_basic_grid_new(char const *name, struct yac_basic_grid_data grid_data)
struct yac_basic_grid * yac_basic_grid_empty_new(char const *name)
void yac_basic_grid_delete(struct yac_basic_grid *grid)
size_t yac_basic_grid_add_mask_nocpy(struct yac_basic_grid *grid, enum yac_location location, int const *mask, char const *mask_name)
enum yac_interp_ncc_weight_type weight_type
int main(int argc, char **argv)
void yac_dist_grid_pair_delete(struct yac_dist_grid_pair *grid_pair)
struct yac_dist_grid_pair * yac_dist_grid_pair_new(struct yac_basic_grid *grid_a, struct yac_basic_grid *grid_b, MPI_Comm comm)
void yac_interp_grid_delete(struct yac_interp_grid *interp_grid)
struct yac_interp_grid * yac_interp_grid_new(struct yac_dist_grid_pair *grid_pair, char const *src_grid_name, char const *tgt_grid_name, size_t num_src_fields, struct yac_interp_field const *src_fields, struct yac_interp_field const tgt_field)
void yac_interp_method_delete(struct interp_method **method)
struct yac_interp_weights * yac_interp_method_do_search(struct interp_method **method, struct yac_interp_grid *interp_grid)
struct interp_method * yac_interp_method_file_new(char const *weight_file_name)
void yac_interp_weights_delete(struct yac_interp_weights *weights)
void yac_interp_weights_write_to_file(struct yac_interp_weights *weights, char const *filename, char const *src_grid_name, char const *tgt_grid_name, size_t src_grid_size, size_t tgt_grid_size)
struct yac_interpolation * yac_interp_weights_get_interpolation(struct yac_interp_weights *weights, enum yac_interp_weights_reorder_type reorder, size_t collection_size, double frac_mask_fallback_value, double scaling_factor, double scaling_summand, char const *yaxt_exchanger_name)
yac_interp_weights_reorder_type
@ YAC_MAPPING_ON_TGT
weights will be applied at target processes
@ YAC_MAPPING_ON_SRC
weights will be appied at source processes
void yac_interpolation_delete(struct yac_interpolation *interp)
void yac_interpolation_execute_get(struct yac_interpolation *interp, double **tgt_field)
void yac_interpolation_execute_put(struct yac_interpolation *interp, double ***src_fields)
double const YAC_FRAC_MASK_NO_VALUE
enum yac_location location
struct yac_interp_field tgt_field
struct yac_dist_grid_pair * grid_pair
struct yac_interp_field src_fields[]