test_dummy_coupling10_c.c

Go to the documentation of this file.

// Copyright (c) 2024 The YAC Authors
//
// SPDX-License-Identifier: BSD-3-Clause
 
#include <stdlib.h>
#include <stdio.h>
#include <mpi.h>
 
#include <yac.h>
#include "tests.h"
#include "geometry.h"
 
int main(void) {
 
  yac_cinit();
  yac_cdef_calendar(YAC_PROLEPTIC_GREGORIAN);
  yac_cdef_datetime("1990-01-01T00:00:00", "2000-01-01T00:00:00");
 
  int comm_rank, comm_size;
  MPI_Comm_rank(MPI_COMM_WORLD, &comm_rank);
  MPI_Comm_size(MPI_COMM_WORLD, &comm_size);
 
  enum {NUM_PROCS = 4};
  if (comm_size != NUM_PROCS) {
    if (comm_rank == 0) PUT_ERR("This test must be run with 4 processes");
    yac_cfinalize();
    return TEST_EXIT_CODE;
  }
 
  // Define component (setup only has a single component)
  char const * comp_name = "main_component";
  int comp_id;
  yac_cdef_comp(comp_name, &comp_id);
 
  // Global grid dimensions
  enum {
    GLOBAL_NUM_CELLS_X = 16,
    GLOBAL_NUM_CELLS_Y = 16,
  };
 
  // Allocate and fill global coordinates for grid
  double vertices_x[GLOBAL_NUM_CELLS_X+1];
  double vertices_y[GLOBAL_NUM_CELLS_Y+1];
  double cells_x[GLOBAL_NUM_CELLS_X];
  double cells_y[GLOBAL_NUM_CELLS_Y];
  for (size_t i = 0; i <= GLOBAL_NUM_CELLS_X; ++i)
    vertices_x[i] =
      (-5.0 + (10.0 * (double)i) / (double)GLOBAL_NUM_CELLS_X) * YAC_RAD;
  for (size_t j = 0; j <= GLOBAL_NUM_CELLS_Y; ++j)
    vertices_y[j] =
      (-5.0 + (10.0 * (double)j) / (double)GLOBAL_NUM_CELLS_Y) * YAC_RAD;
  for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i)
    cells_x[i] = 0.5 * (vertices_x[i] + vertices_x[i+1]);
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j)
    cells_y[j] = 0.5 * (vertices_y[j] + vertices_y[j+1]);
 
  // Define local grid regular grid (each process uses one fourth of grid rows)
  int local_start_x = 0;
  int local_start_y = (comm_rank * GLOBAL_NUM_CELLS_Y) / comm_size;
  int local_count_x = GLOBAL_NUM_CELLS_X;
  int local_count_y =
    ((comm_rank + 1) * GLOBAL_NUM_CELLS_Y) / comm_size - local_start_y;
  char const * grid_name = "main_grid";
  int nbr_cells[2] = {local_count_x, local_count_y};
  int nbr_vertices[2] = {local_count_x + 1, local_count_y + 1};
  int cyclic[2] = {0,0};
  int grid_id;
  yac_cdef_grid_reg2d(
    grid_name, nbr_vertices, cyclic,
    vertices_x + local_start_x, vertices_y + local_start_y, &grid_id);
 
  // Define cell center coordinates
  int cell_point_id;
  yac_cdef_points_reg2d(
    grid_id, nbr_cells, YAC_LOCATION_CELL,
    cells_x + local_start_x, cells_y + local_start_y, &cell_point_id);
 
  // Generate global coast and ocean mask
  enum {COAST_COLUMN = 8}; // column 8 is the coast
  int coast_mask[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
  int ocean_mask[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j) {
    for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i) {
 
      // column 8 is the coast
      coast_mask[j][i] = (i == COAST_COLUMN);
 
      // everything right of coast is ocean
      ocean_mask[j][i] = (i > COAST_COLUMN);
    }
  }
  int coast_mask_id, ocean_mask_id;
  yac_cdef_mask(
    grid_id, local_count_x * local_count_y, YAC_LOCATION_CELL,
    &coast_mask[local_start_y][0], &coast_mask_id);
  yac_cdef_mask(
    grid_id, local_count_x * local_count_y, YAC_LOCATION_CELL,
    &ocean_mask[local_start_y][0], &ocean_mask_id);
 
  // Define fields
  enum {COLLECTION_SIZE = 1, NUM_POINTSETS = 1};
  char const * timestep = "1";
  int temp_field_id, sst_field_id, runoff_coast_field_id, runoff_ocean_field_id;
  // temp: no mask
  yac_cdef_field(
    "temp", comp_id, &cell_point_id, NUM_POINTSETS,
    COLLECTION_SIZE, timestep, YAC_TIME_UNIT_HOUR, &temp_field_id);
  // sst: ocean mask
  yac_cdef_field_mask(
    "sst", comp_id, &cell_point_id, &ocean_mask_id, NUM_POINTSETS,
    COLLECTION_SIZE, timestep, YAC_TIME_UNIT_HOUR, &sst_field_id);
  // runoff: coast mask
  yac_cdef_field_mask(
    "runoff", comp_id, &cell_point_id, &coast_mask_id, NUM_POINTSETS,
    COLLECTION_SIZE, timestep, YAC_TIME_UNIT_HOUR, &runoff_coast_field_id);
  // runoff: ocean mask
  yac_cdef_field_mask(
    "river_runoff", comp_id, &cell_point_id, &ocean_mask_id, NUM_POINTSETS,
    COLLECTION_SIZE, timestep, YAC_TIME_UNIT_HOUR, &runoff_ocean_field_id);
 
 
  // Define coupling
  // temp -> sst using 1NN
  enum {SRC_LAG = 0, TGT_LAG = 0};
  int interp_stack_1nn;
  yac_cget_interp_stack_config_from_string_json(
    "[{\"nnn\": {\"n\": 1}}]", &interp_stack_1nn);
  yac_cdef_couple(
    comp_name, grid_name, "temp",
    comp_name, grid_name, "sst",
    timestep, YAC_TIME_UNIT_HOUR, YAC_REDUCTION_TIME_NONE,
    interp_stack_1nn, SRC_LAG, TGT_LAG);
  yac_cfree_interp_stack_config(interp_stack_1nn);
  // runoff (coast) -> runoff (ocean) using SPMAP
  int interp_stack_spmap;
  yac_cget_interp_stack_config_from_string_json(
    "[\"source_to_target_map\", {\"fixed\": {\"user_value\": 999}}]",
    &interp_stack_spmap);
  yac_cdef_couple(
    comp_name, grid_name, "runoff",
    comp_name, grid_name, "river_runoff",
    timestep, YAC_TIME_UNIT_HOUR, YAC_REDUCTION_TIME_NONE,
    interp_stack_spmap, SRC_LAG, TGT_LAG);
  yac_cfree_interp_stack_config(interp_stack_spmap);
 
  // End of definition
  yac_cenddef();
 
 
  // Set up dummy data for fields
  double global_temp_data[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
  double global_runoff_data[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
 
  // Initialize temp: simple pattern, e.g., temp = 100.0 + 10.0 * j + i
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j) {
    for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i) {
      global_temp_data[j][i] = 100.0 + 10.0 * (double)j + (double)i;
    }
  }
 
  // Initialize runoff: set coast column to 42.0, rest to -1.0
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j) {
    for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i) {
      if (i == COAST_COLUMN) {
        global_runoff_data[j][i] = (double)j;
      } else {
        global_runoff_data[j][i] = -1.0;
      }
    }
  }
 
  // Output buffers for received fields
  double sst_data[local_count_y][local_count_x];
  double river_runoff_data[local_count_y][local_count_x];
  for (int j = 0; j < local_count_y; ++j) {
    for (int i = 0; i < local_count_x; ++i) {
      sst_data[j][i] = -2.0;
      river_runoff_data[j][i] = -2.0;
    }
  }
 
  // Do single exchange of fields using yac_cexchange (pairwise)
  int info, ierror = 0;
  yac_cexchange_(
    temp_field_id, sst_field_id, COLLECTION_SIZE,
    &global_temp_data[local_start_y][0], &sst_data[0][0],
    &info, &info, &ierror);
  yac_cexchange_(
    runoff_coast_field_id, runoff_ocean_field_id, COLLECTION_SIZE,
    &global_runoff_data[local_start_y][0], &river_runoff_data[0][0],
    &info, &info, &ierror);
 
  // Reference for sst:
  // for ocean points, should be nearest temp (column > COAST_COLUMN)
  double ref_global_sst_data[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j) {
    for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i) {
      if (i > COAST_COLUMN) {
        // 1NN: nearest is itself (since grid is regular)
        ref_global_sst_data[j][i] = global_temp_data[j][i];
      } else {
        ref_global_sst_data[j][i] = -2.0; // not ocean
      }
    }
  }
 
  // Reference for river_runoff:
  // ocean points adjacent to coast should receive its value,
  // other ocean points the fixed value
  double ref_global_river_runoff_data[GLOBAL_NUM_CELLS_Y][GLOBAL_NUM_CELLS_X];
  for (size_t j = 0; j < GLOBAL_NUM_CELLS_Y; ++j) {
    for (size_t i = 0; i < GLOBAL_NUM_CELLS_X; ++i) {
      if (i == COAST_COLUMN + 1) {
        // SPMAP: ocean receives from coast column in same row
        ref_global_river_runoff_data[j][i] = global_runoff_data[j][COAST_COLUMN];
      } else if (i > COAST_COLUMN + 1) {
        ref_global_river_runoff_data[j][i] = 999.0; // fixed value
      } else {
        ref_global_river_runoff_data[j][i] = -2.0; // not ocean
      }
    }
  }
 
  // Check received sst and river_runoff
  for (int j = local_start_y; j < local_count_y; ++j) {
    for (int i = local_start_x; i < local_count_x; ++i) {
      if (sst_data[j][i] !=
          ref_global_sst_data[local_start_y + j][local_start_x + i]) {
        PUT_ERR("Error in sst data");
      }
      if (river_runoff_data[j][i] !=
          ref_global_river_runoff_data[local_start_y + j][local_start_x + i]) {
        PUT_ERR("Error in sst data");
      }
    }
  }
 
  // Finalize YAC
  yac_cfinalize ();
 
  return TEST_EXIT_CODE;
}