Skip to content

Commit

Permalink
fix: fix calculation of source index for winds using yac edges fields
Browse files Browse the repository at this point in the history
  • Loading branch information
yoctoyotta1024 committed Sep 9, 2024
1 parent 41e47fb commit 4fad417
Show file tree
Hide file tree
Showing 2 changed files with 91 additions and 34 deletions.
108 changes: 79 additions & 29 deletions libs/coupldyn_yac/yac_cartesian_dynamics.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
* Author: Wilton Loch (WL)
* Additional Contributors: Clara Bayley (CB)
* -----
* Last Modified: Saturday 24th August 2024
* Last Modified: Monday 9th September 2024
* Modified By: CB
* -----
* License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -132,48 +132,98 @@ void create_grid_and_points_definitions(const Config &config, const std::array<s
}

/*
fill's target_array with values from yac_raw_data at multiplied by their conversion factor
fill's target_array with values from yac_raw_data at multiplied by their
conversion factor. Special handling when grid_points = 1 or 2 such that only
edge grid points for every other layer are used (to correctly choose only edges
in that direction and not another when using YAC's edge point definition).
*/
void CartesianDynamics::receive_yac_field(unsigned int yac_field_id, double **yac_raw_data,
void CartesianDynamics::receive_yac_field(const unsigned int grid_points,
const unsigned int yac_field_id, double **yac_raw_data,
std::vector<double> &target_array,
const size_t ndims_north, const size_t ndims_east,
const size_t vertical_levels,
double conversion_factor = 1.0) const {
const double conversion_factor = 1.0) const {
int info, error;
yac_cget(yac_field_id, vertical_levels, yac_raw_data, &info, &error);

for (size_t j = 0; j < ndims_north; j++) {
for (size_t i = 0; i < ndims_east; i++) {
for (size_t k = 0; k < vertical_levels; k++) {
auto vertical_idx = k;
auto source_idx = j * ndims_east + i;
auto ii = (ndims_east * j + i) * vertical_levels + k;
target_array[ii] = yac_raw_data[vertical_idx][source_idx] / conversion_factor;
switch (grid_points) {
case 0:
/* handles fields defined on centres of cells of horizontal
(2-D) grid (grid_points = CENTRES), e.g. temp or wvel
*/
for (size_t j = 0; j < ndims_north; j++) {
for (size_t i = 0; i < ndims_east; i++) {
for (size_t k = 0; k < vertical_levels; k++) {
auto vertical_idx = k;
auto source_idx = ndims_east * j + i;
auto ii = (ndims_east * j + i) * vertical_levels + k;
target_array[ii] = yac_raw_data[vertical_idx][source_idx] / conversion_factor;
}
}
}
}
return;
case 1:
/* handles fields defined on longitude edges of cells of
horizontal grid (grid_points = LONGITUDE_EDGES), e.g uvel
*/
for (size_t j = 0; j < ndims_north; j++) {
for (size_t i = 0; i < ndims_east; i++) {
for (size_t k = 0; k < vertical_levels; k++) {
auto vertical_idx = k;
auto edges_ndims_east = ndims_east * 2 - 1;
auto source_idx = edges_ndims_east * j + ndims_east - 1 + i;
auto ii = (ndims_east * j + i) * vertical_levels + k;
target_array[ii] = yac_raw_data[vertical_idx][source_idx] / conversion_factor;
}
}
}
return;

case 2:
/* handles fields defined on latitude edges of cells of
horizontal grid (grid_points = LATITUDE_EDGES), e.g. vvel
*/
for (size_t j = 0; j < ndims_north; j++) {
for (size_t i = 0; i < ndims_east; i++) {
for (size_t k = 0; k < vertical_levels; k++) {
auto vertical_idx = k;
auto edges_ndims_east = ndims_east * 2 + 1;
auto source_idx = edges_ndims_east * j + i;
auto ii = (ndims_east * j + i) * vertical_levels + k;
target_array[ii] = yac_raw_data[vertical_idx][source_idx] / conversion_factor;
}
}
}
return;
}
}

/* This subroutine is the main entry point for receiving data from YAC.
* It checks the dimensionality of the simulation based on the config data. */
void CartesianDynamics::receive_fields_from_yac() {
receive_yac_field(temp_yac_id, yac_raw_cell_data, temp, ndims[NORTHWARD], ndims[EASTWARD],
ndims[VERTICAL], dlc::TEMP0);
receive_yac_field(pressure_yac_id, yac_raw_cell_data, press, ndims[NORTHWARD], ndims[EASTWARD],
ndims[VERTICAL], dlc::P0);
receive_yac_field(qvap_yac_id, yac_raw_cell_data, qvap, ndims[NORTHWARD], ndims[EASTWARD],
ndims[VERTICAL]);
receive_yac_field(qcond_yac_id, yac_raw_cell_data, qcond, ndims[NORTHWARD], ndims[EASTWARD],
ndims[VERTICAL]);

receive_yac_field(vertical_wind_yac_id, yac_raw_vertical_wind_data, wvel, ndims[NORTHWARD],
ndims[EASTWARD], ndims[VERTICAL] + 1, dlc::W0);

receive_yac_field(eastward_wind_yac_id, yac_raw_edge_data, uvel, ndims[NORTHWARD],
ndims[EASTWARD] + 1, ndims[VERTICAL], dlc::W0);

receive_yac_field(northward_wind_yac_id, yac_raw_edge_data, vvel, ndims[NORTHWARD] + 1,
ndims[EASTWARD], ndims[VERTICAL], dlc::W0);
enum grid_points {
CENTRES,
LONGITUDE_EDGES,
LATITUDE_EDGES,
};

receive_yac_field(CENTRES, temp_yac_id, yac_raw_cell_data, temp, ndims[NORTHWARD],
ndims[EASTWARD], ndims[VERTICAL], dlc::TEMP0);
receive_yac_field(CENTRES, pressure_yac_id, yac_raw_cell_data, press, ndims[NORTHWARD],
ndims[EASTWARD], ndims[VERTICAL], dlc::P0);
receive_yac_field(CENTRES, qvap_yac_id, yac_raw_cell_data, qvap, ndims[NORTHWARD],
ndims[EASTWARD], ndims[VERTICAL]);
receive_yac_field(CENTRES, qcond_yac_id, yac_raw_cell_data, qcond, ndims[NORTHWARD],
ndims[EASTWARD], ndims[VERTICAL]);

receive_yac_field(CENTRES, vertical_wind_yac_id, yac_raw_vertical_wind_data, wvel,
ndims[NORTHWARD], ndims[EASTWARD], ndims[VERTICAL] + 1, dlc::W0);

receive_yac_field(LONGITUDE_EDGES, eastward_wind_yac_id, yac_raw_edge_data, uvel,
ndims[NORTHWARD], ndims[EASTWARD] + 1, ndims[VERTICAL], dlc::W0);

receive_yac_field(LATITUDE_EDGES, northward_wind_yac_id, yac_raw_edge_data, vvel,
ndims[NORTHWARD] + 1, ndims[EASTWARD], ndims[VERTICAL], dlc::W0);
}

CartesianDynamics::CartesianDynamics(const Config &config, const std::array<size_t, 3> i_ndims,
Expand Down
17 changes: 12 additions & 5 deletions libs/coupldyn_yac/yac_cartesian_dynamics.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
* Author: Clara Bayley (CB)
* Additional Contributors:
* -----
* Last Modified: Saturday 24th August 2024
* Last Modified: Wednesday 4th September 2024
* Modified By: CB
* -----
* License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -123,10 +123,17 @@ struct CartesianDynamics {
/* Public call to receive data from YAC
* If the problem is 2D turns into a wrapper for receive_hor_slice_from_yac */
void receive_fields_from_yac();
void receive_yac_field(unsigned int yac_field_id, double **yac_raw_data,
std::vector<double> &target_array, const size_t ndims_north,
const size_t ndims_east, const size_t vertical_levels,
double conversion_factor) const;

/*
fill's target_array with values from yac_raw_data at multiplied by their
conversion factor. Special handling when grid_points = 1 or 2 such that only
edge grid points for every other layer are used (to correctly choose only edges
in that direction and not another when using YAC's edge point definition).
*/
void receive_yac_field(const unsigned int grid_points, const unsigned int yac_field_id,
double **yac_raw_data, std::vector<double> &target_array,
const size_t ndims_north, const size_t ndims_east,
const size_t vertical_levels, const double conversion_factor) const;
};

/* type satisfying CoupledDyanmics solver concept
Expand Down

0 comments on commit 4fad417

Please sign in to comment.