Commit d5e42ee4 authored by Moritz Hanke's avatar Moritz Hanke

improves performance of local search in spmap interpolation method

parent 6005eb47
......@@ -32,12 +32,14 @@
*/
#include <string.h>
#include <math.h>
#include <sys/time.h>
#include "interpolation_method.h"
#include "interpolation_method_spmap.h"
#include "ensure_array_size.h"
#include "geometry.h"
#include "grid.h"
#include "sphere_part.h"
static void
do_search_interpolation_method_spmap(struct interpolation_method ** method,
......@@ -135,6 +137,8 @@ struct search_data {
struct src_point_search_data * src_point_data;
unsigned num_src_points;
double * x_coordinates;
double * y_coordinates;
};
struct tgt_point_data {
......@@ -142,6 +146,7 @@ struct tgt_point_data {
double * y_coordinates;
int * global_ids;
unsigned num_tgt_points;
struct point_sphere_part_search * point_sphere_part;
};
struct halo_cells {
......@@ -486,6 +491,18 @@ init_local_search_data(struct search_data * search_data, struct domain * domain,
search_data->src_point_data =
realloc(src_point_data, num_src_points * sizeof(*src_point_data));
search_data->num_src_points = num_src_points;
search_data->x_coordinates =
malloc(num_src_points * sizeof(*(search_data->x_coordinates)));
search_data->y_coordinates =
malloc(num_src_points * sizeof(*(search_data->y_coordinates)));
for (unsigned i = 0; i < num_src_points; ++i) {
search_data->x_coordinates[i] =
search_data->src_point_data[i].circle.base_point[0];
search_data->y_coordinates[i] =
search_data->src_point_data[i].circle.base_point[1];
}
}
// searches through local target points
......@@ -496,34 +513,58 @@ static void search_local(struct search_data * search_data,
struct src_point_search_data * src_point_data = search_data->src_point_data;
unsigned num_src_points = search_data->num_src_points;
for (unsigned tgt_idx = 0; tgt_idx < tgt_points.num_tgt_points; ++tgt_idx) {
unsigned * local_point_ids = NULL;
unsigned local_point_ids_array_size = 0;
unsigned * num_local_point_ids =
malloc(search_data->num_src_points * sizeof(*num_local_point_ids));
yac_point_sphere_part_search_NN(tgt_points.point_sphere_part,
search_data->num_src_points,
search_data->x_coordinates,
search_data->y_coordinates,
&local_point_ids, &local_point_ids_array_size,
num_local_point_ids);
// if more than on local point, use highest; because tgt point are sorted by
// global_id, it will also be the highest one
unsigned offset = 0;
for (unsigned src_idx = 0; src_idx < search_data->num_src_points; ++src_idx) {
if (num_local_point_ids[src_idx] == 0) continue;
unsigned tgt_idx = local_point_ids[offset];
double tgt_coords_xyz[3];
LLtoXYZ(tgt_points.x_coordinates[tgt_idx],
tgt_points.y_coordinates[tgt_idx], tgt_coords_xyz);
double distance =
get_vector_angle(tgt_coords_xyz,
src_point_data[src_idx].circle.base_vector);
// if the result is not better than the one already found
if (distance > src_point_data[src_idx].circle.inc_angle) continue;
int tgt_global_id = tgt_points.global_ids[tgt_idx];
// find the target point with the highest global id
for (unsigned i = 1; i < num_local_point_ids[src_idx]; ++i) {
unsigned tmp_tgt_idx = local_point_ids[offset + i];
int tmp_tgt_global_id = tgt_points.global_ids[tmp_tgt_idx];
if (tmp_tgt_global_id > tgt_global_id) {
tgt_idx = tmp_tgt_idx;
tgt_global_id = tmp_tgt_global_id;
}
}
for (unsigned src_idx = 0; src_idx < search_data->num_src_points;
++src_idx) {
double distance =
get_vector_angle(tgt_coords_xyz,
src_point_data[src_idx].circle.base_vector);
// if the distance is bigger than the one already found
if (distance > src_point_data[src_idx].circle.inc_angle) continue;
// if the distance for the current target point is the same as the one
// already found and the global index of the target points is smaller
if ((distance == src_point_data[src_idx].circle.inc_angle) &&
(src_point_data[src_idx].tgt_point.global_id >=
tgt_points.global_ids[tgt_idx])) continue;
src_point_data[src_idx].circle.inc_angle = distance;
src_point_data[src_idx].tgt_point.global_id = tgt_global_id;
src_point_data[src_idx].tgt_point.rank = rank;
src_point_data[src_idx].circle.inc_angle = distance;
src_point_data[src_idx].tgt_point.global_id =
tgt_points.global_ids[tgt_idx];
src_point_data[src_idx].tgt_point.rank = rank;
}
offset += num_local_point_ids[src_idx];
}
free(local_point_ids);
free(num_local_point_ids);
}
// gets the bounding circle and the owner rank of all halo cells
......@@ -969,6 +1010,9 @@ do_search_final(struct intermediate_method_data * data,
free(data->search_data.src_point_data[i].checked_ranks);
free(data->search_data.src_point_data[i].to_be_checked_ranks);
}
yac_delete_point_sphere_part_search(data->tgt_points.point_sphere_part);
free(data->search_data.y_coordinates);
free(data->search_data.x_coordinates);
free(data->search_data.src_point_data);
free(data);
}
......@@ -1763,6 +1807,8 @@ pack_and_send_search_request_answer(
yac_comm_bsend(buffer, buffer_size, COMM_PACKED, remote_rank, data_tag, comm);
free(buffer);
free(search_data->y_coordinates);
free(search_data->x_coordinates);
free(search_data->src_point_data);
}
......@@ -1801,6 +1847,11 @@ recv_search_request_cont_callback(struct communicator * comm,
int position = 0;
search_data.x_coordinates =
malloc(search_data.num_src_points * sizeof(*(search_data.x_coordinates)));
search_data.y_coordinates =
malloc(search_data.num_src_points * sizeof(*(search_data.y_coordinates)));
// unpack and initialise search_data
for (unsigned i = 0; i < search_data.num_src_points; ++i) {
......@@ -1813,6 +1864,10 @@ recv_search_request_cont_callback(struct communicator * comm,
yac_comm_unpack(msg.buffer, msg.count, &position,
search_data.src_point_data[i].circle.base_point, 2,
COMM_DBLE, comm);
search_data.x_coordinates[i] =
search_data.src_point_data[i].circle.base_point[0];
search_data.y_coordinates[i] =
search_data.src_point_data[i].circle.base_point[1];
yac_comm_unpack(msg.buffer, msg.count, &position,
&(search_data.src_point_data[i].circle.inc_angle), 1,
COMM_DBLE, comm);
......@@ -1932,6 +1987,10 @@ do_search_interpolation_method_spmap(
// sort target points by global_ids
sort_target_points(global_ids, x_coordinates, y_coordinates, num_points);
// generate search tree for the target points
data->tgt_points.point_sphere_part =
yac_point_sphere_part_search_new(num_points, x_coordinates, y_coordinates);
// get all non-masked local source core points and sort them by their
// global ids
init_local_search_data(&(data->search_data), method_spmap->domain,
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment