Commit 2a6751f7 authored by Francesco Tomba's avatar Francesco Tomba
Browse files

Fix MPI count overflow, thread-safety and misc warnings

- Use MPI_Type_contiguous derived types for point exchange and kNN
  alltoallv blocks so int counts stay in element units (BUG_REPORT 3.6)
- Move LOG_WRITE-hidden barriers out of NDEBUG; fix TIME_STOP/LOG_WRITE
- Fix RMA window epochs in Heuristic1/2 (remove leading/trailing fences)
- H1 zero-init loop bound mpi_rank -> world_size (BUG_REPORT 3.7)
- big_ordered_buffer_to_file: receive via probed count_recv (3.4)
- Gather size_t counts as MPI_UINT64_T (3.5)
- Replace MPI_Testall misuse with MPI_Waitall (3.3)
- Deduplicate MY_MALLOC; single definition in heap.h
- Thread-safe dadp_rand (rand_r) replacing rand()/srand in tree code
- getopt: required_argument for k/z/knn-strategy, strtol validation
- fopen failure hardening with exit(1)
- Remove dead code: va_start redefinitions, mock read args, data_dims,
  bits/time.h; fix %lu/%u, MPI_UNSIGNED_LONG, void* arithmetic
parent 585e3063
Loading
Loading
Loading
Loading
+21 −20
Original line number Diff line number Diff line
#include "adp.h"
#include "mpi.h"
#include <bits/time.h>
#include <time.h>
#include <stdio.h>
#include <string.h>
#include <time.h>
@@ -166,7 +166,7 @@ void compute_density_kstarnn_rma_v2(global_context_t* ctx, const float_t d, int
    #pragma omp parallel for
    for(idx_t i = 0; i < ctx -> local_n_points; ++i)
    {
        for(idx_t k = 0; k <= ctx -> k; ++k)
        for(idx_t k = 0; k < ctx -> k; ++k)
        {
            local_datapoints[i].ngbh[k].value = omega * pow(local_datapoints[i].ngbh[k].value, d/2.);  
        }
@@ -702,8 +702,8 @@ clusters_t Heuristic1(global_context_t *ctx)
    MPI_Win win_datapoints;
    MPI_Win_create(ctx -> local_datapoints, ctx -> local_n_points * sizeof(datapoint_info_t), 
                   1, MPI_INFO_NULL, ctx -> mpi_communicator, &win_datapoints);
    MPI_Win_fence(0, win_datapoints);
    MPI_Win_lock_all(0,  win_datapoints);
    MPI_Barrier(ctx -> mpi_communicator);
    //
#if !defined(THREAD_FUNNELED)
    #pragma omp parallel for
@@ -783,6 +783,8 @@ clusters_t Heuristic1(global_context_t *ctx)

    elapsed_time = TIME_STOP;
    LOG_WRITE("Putative centers", elapsed_time);
    MPI_Win_sync(win_datapoints);
    MPI_Barrier(ctx -> mpi_communicator);

    TIME_START;

@@ -965,6 +967,7 @@ clusters_t Heuristic1(global_context_t *ctx)
    elapsed_time = TIME_STOP;
    LOG_WRITE("Finding centers to prune", elapsed_time);
    TIME_START;
    MPI_Barrier(ctx -> mpi_communicator);
    
    idx_t tot_removal = 0;
    for(idx_t p = 0; p < n; ++p)
@@ -1005,7 +1008,7 @@ clusters_t Heuristic1(global_context_t *ctx)

    //zero_out_all

    for(int i = 0; i < ctx -> mpi_rank; ++i)
    for(int i = 0; i < ctx -> world_size; ++i)
    {
        recv_counts[i] = 0;
        send_counts[i] = 0;
@@ -1095,13 +1098,11 @@ clusters_t Heuristic1(global_context_t *ctx)
    */


    for(int i = 0; i < ctx -> world_size; ++i)
    {
        recv_counts[i] = recv_counts[i] * sizeof(center_removal_t);
        send_counts[i] = send_counts[i] * sizeof(center_removal_t);
        recv_displs[i] = recv_displs[i] * sizeof(center_removal_t);
        send_displs[i] = send_displs[i] * sizeof(center_removal_t);
    }
    /* keep counts in units of center_removal_t via a derived type so that
     * count * sizeof(center_removal_t) never overflows the int counts */
    MPI_Datatype mpi_center_removal_type;
    MPI_Type_contiguous(sizeof(center_removal_t), MPI_BYTE, &mpi_center_removal_type);
    MPI_Type_commit(&mpi_center_removal_type);

    //allocate buffer to recieve center elminiations
    
@@ -1109,14 +1110,17 @@ clusters_t Heuristic1(global_context_t *ctx)
    center_removal_t* recv_removals = (center_removal_t*)MY_MALLOC(tot_recv_counts * sizeof(center_removal_t));

    // all to all
    MPI_Alltoallv(removal_buffer, send_counts, send_displs, MPI_BYTE, 
                   recv_removals, recv_counts, recv_displs, MPI_BYTE, ctx -> mpi_communicator);
    MPI_Alltoallv(removal_buffer, send_counts, send_displs, mpi_center_removal_type, 
                   recv_removals, recv_counts, recv_displs, mpi_center_removal_type, ctx -> mpi_communicator);

    MPI_Type_free(&mpi_center_removal_type);

    // merge into the mask

    elapsed_time = TIME_STOP;
    LOG_WRITE("Communicating eliminations", elapsed_time);
    TIME_START;
    MPI_Barrier(ctx -> mpi_communicator);
    
    #pragma omp parallel for schedule(dynamic)
    for(idx_t i = 0; i < tot_recv_counts; ++i)
@@ -1203,6 +1207,7 @@ clusters_t Heuristic1(global_context_t *ctx)
    elapsed_time = TIME_STOP;
    LOG_WRITE("Merging", elapsed_time);

    MPI_Barrier(ctx -> mpi_communicator);
    TIME_START;

    int n_centers = (int)actual_centers.count;
@@ -1348,10 +1353,8 @@ clusters_t Heuristic1(global_context_t *ctx)
    }

    MPI_Win_unlock_all(win_datapoints);
    MPI_Win_fence(0, win_datapoints);
    MPI_Win_free(&win_datapoints);

    MPI_Barrier(ctx -> mpi_communicator);
    MPI_Win_free(&win_datapoints);

    free(dp_info_ptrs);
    free(max_rho.data);
@@ -1370,6 +1373,7 @@ clusters_t Heuristic1(global_context_t *ctx)

    elapsed_time = TIME_STOP;
    LOG_WRITE("Cluster assign", elapsed_time);
    MPI_Barrier(ctx -> mpi_communicator);

    free(actual_centers.data);
    actual_centers.size  = tot_centers;
@@ -1407,9 +1411,6 @@ void Heuristic2(global_context_t* ctx, clusters_t* cluster)
    MPI_Win_create(ctx -> local_datapoints, ctx -> local_n_points * sizeof(datapoint_info_t), 1, MPI_INFO_NULL, ctx -> mpi_communicator, &dp_info_win);
    MPI_Win_create(ctx -> __local_heap_buffers, ctx -> local_n_points * ctx -> k * sizeof(heap_node_t), 1, MPI_INFO_NULL, ctx -> mpi_communicator, &ngbh_win);

    MPI_Win_fence(0, dp_info_win);
    MPI_Win_fence(0, ngbh_win);

    MPI_Win_lock_all(0, dp_info_win);
    MPI_Win_lock_all(0, ngbh_win);
    MPI_Barrier(ctx -> mpi_communicator);
+24 −16
Original line number Diff line number Diff line
@@ -4,7 +4,7 @@
#include <time.h>
#include "../adp/adp.h"

#define ARRAY_INCREMENT 100
#define LU_ARRAY_INCREMENT 100

#define FREE_NOT_NULL(x) if(x){free(x); x = NULL;}
#define PATH_LEN 500
@@ -241,9 +241,9 @@ void generate_random_matrix(

void lu_dynamic_array_allocate(lu_dynamic_array_t * a)
{
    a -> data = (idx_t*)malloc(ARRAY_INCREMENT*sizeof(idx_t));
    a -> data = (idx_t*)malloc(LU_ARRAY_INCREMENT*sizeof(idx_t));
    a -> count = 0;
    a -> size = ARRAY_INCREMENT;
    a -> size = LU_ARRAY_INCREMENT;
}

void lu_dynamic_array_pushBack(lu_dynamic_array_t * a, idx_t p)
@@ -254,7 +254,7 @@ void lu_dynamic_array_pushBack(lu_dynamic_array_t * a, idx_t p)
        a -> count += 1;
    }
    else{
        a -> size += ARRAY_INCREMENT;
        a -> size += LU_ARRAY_INCREMENT;
        a -> data = realloc(a -> data, a -> size * sizeof(idx_t));
        a -> data[a -> count] =  p;
        a -> count += 1;
@@ -332,7 +332,7 @@ float_t* read_data_file(global_context_t *ctx, const char *fname, const idx_t nd
        ctx -> dims = ndims;
        ctx -> n_points = n / ctx -> dims;

        mpi_printf(ctx, "Got ndims %lu npoints %lu\n", ctx -> dims, ctx -> n_points);
        mpi_printf(ctx, "Got ndims %u npoints %lu\n", (unsigned)ctx -> dims, ctx -> n_points);
        fclose(f);

        for (uint64_t i = 0; i < n; ++i) data[i] = (float_t)(df[i]);
@@ -351,7 +351,7 @@ float_t* read_data_file(global_context_t *ctx, const char *fname, const idx_t nd
        ctx -> dims = ndims;
        ctx -> n_points = n / ctx -> dims;

        mpi_printf(ctx, "Got ndims %lu npoints %lu\n", ctx -> dims, ctx -> n_points);
        mpi_printf(ctx, "Got ndims %u npoints %lu\n", (unsigned)ctx -> dims, ctx -> n_points);
        fclose(f);

        for (uint64_t i = 0; i < n; ++i) data[i] = (float_t)(df[i]);
@@ -508,22 +508,24 @@ void ordered_buffer_to_file(global_context_t* ctx, void* buffer, size_t el_size,
    uint64_t tot_n = 0;
    MPI_Reduce(&n, &tot_n, 1, MPI_UINT64_T , MPI_SUM, 0, ctx -> mpi_communicator);

    uint64_t* ppp64 = NULL;
    if(I_AM_MASTER) 
    {
        tmp_data = (void*)MY_MALLOC(el_size * tot_n );
        ppp      = (int*)MY_MALLOC(ctx -> world_size * sizeof(int));
        displs   = (int*)MY_MALLOC(ctx -> world_size * sizeof(int));
        ppp64    = (uint64_t*)MY_MALLOC(ctx -> world_size * sizeof(uint64_t));

    }
    
    int nn = (int)n;
    MPI_Gather(&nn, 1, MPI_INT, ppp, 1, MPI_INT, 0, ctx -> mpi_communicator);
    MPI_Gather(&n, 1, MPI_UINT64_T, ppp64, 1, MPI_UINT64_T, 0, ctx -> mpi_communicator);

    if(I_AM_MASTER)
    {
        displs[0] = 0;
        for(int i = 0; i < ctx -> world_size; ++i) ppp[i]    = el_size  * ppp[i];
        for(int i = 0; i < ctx -> world_size; ++i) ppp[i]    = (int)(el_size  * ppp64[i]);
        for(int i = 1; i < ctx -> world_size; ++i) displs[i] = displs[i - 1] + ppp[i - 1];
        free(ppp64);
            
    }

@@ -536,6 +538,7 @@ void ordered_buffer_to_file(global_context_t* ctx, void* buffer, size_t el_size,
        if(!file)
        {
            printf("Cannot open file %s ! Aborting \n", fname);
            exit(1);
        }
        fwrite(tmp_data, 1, el_size * tot_n, file);
        fclose(file);
@@ -556,12 +559,13 @@ void distributed_buffer_to_file(global_context_t* ctx, void* buffer, size_t el_s
    if(!file)
    {
        printf("Cannot open file %s ! Aborting \n", fname);
        exit(1);
    }
    else
    {
        fwrite(buffer, 1, el_size * n, file);
    }
        fclose(file);
    }

}

@@ -617,12 +621,10 @@ void big_ordered_buffer_to_file(global_context_t* ctx, void* buffer, size_t el_s
                MPI_Status status;
                MPI_Probe(r, MPI_ANY_TAG, ctx -> mpi_communicator, &status);

                MPI_Request request;
                int count_recv; 
                int source = status.MPI_SOURCE;
                MPI_Get_count(&status, MPI_BYTE, &count_recv);

                MPI_Recv(tmp_data + displs[r] + already_recv[r], ppp[r], MPI_BYTE, r, r, ctx -> mpi_communicator, MPI_STATUS_IGNORE);
                MPI_Recv((char*)tmp_data + displs[r] + already_recv[r], count_recv, MPI_BYTE, r, r, ctx -> mpi_communicator, MPI_STATUS_IGNORE);
                already_recv[r] += count_recv;

                #ifdef PRINT_ORDERED_BUFFER
@@ -650,6 +652,7 @@ void big_ordered_buffer_to_file(global_context_t* ctx, void* buffer, size_t el_s
        if(!file)
        {
            printf("Cannot open file %s ! Aborting \n", fname);
            exit(1);
        }
        fwrite(tmp_data, 1, el_size * tot_n, file);
        fclose(file);
@@ -670,21 +673,24 @@ void ordered_data_to_file(global_context_t* ctx, const char* fname)
    int* displs;

    MPI_Barrier(ctx -> mpi_communicator);
    uint64_t* ppp64 = NULL;
    if(I_AM_MASTER) 
    {
        tmp_data = (float_t*)MY_MALLOC(ctx -> dims * ctx -> n_points * sizeof(float_t));
        ppp      = (int*)MY_MALLOC(ctx -> world_size * sizeof(int));
        displs   = (int*)MY_MALLOC(ctx -> world_size * sizeof(int));
        ppp64    = (uint64_t*)MY_MALLOC(ctx -> world_size * sizeof(uint64_t));

    }
    
    MPI_Gather(&(ctx -> local_n_points), 1, MPI_INT, ppp, 1, MPI_INT, 0, ctx -> mpi_communicator);
    MPI_Gather(&(ctx -> local_n_points), 1, MPI_UINT64_T, ppp64, 1, MPI_UINT64_T, 0, ctx -> mpi_communicator);

    if(I_AM_MASTER)
    {
        displs[0] = 0;
        for(int i = 0; i < ctx -> world_size; ++i) ppp[i]    = ctx -> dims * ppp[i];
        for(int i = 0; i < ctx -> world_size; ++i) ppp[i]    = (int)(ctx -> dims * ppp64[i]);
        for(int i = 1; i < ctx -> world_size; ++i) displs[i] = displs[i - 1] + ppp[i - 1];
        free(ppp64);
            
    }
    MPI_Gatherv(ctx -> local_data, ctx -> dims * ctx -> local_n_points, 
@@ -701,6 +707,7 @@ void ordered_data_to_file(global_context_t* ctx, const char* fname)
        else
        {
            printf("Cannot open file %s\n", fname);
            exit(1);
        }
        free(tmp_data);
        free(ppp);
@@ -734,13 +741,14 @@ void test_distributed_file_path(global_context_t* ctx, const char* fname)
    if(!file)
    {
        printf("Cannot open file %s ! Aborting \n", out_path_w_proc_name);
        exit(1);
    }
    else
    {
        fprintf(file, "This is only to test if I can open a file in the desidered path\n");
        fprintf(file, "Here will be written the output of dadp\n");
    }
        fclose(file);
    }

}

+6 −23
Original line number Diff line number Diff line
@@ -4,11 +4,11 @@
#include <stdio.h>
#include <mpi.h>
#include <stdint.h>
#include <stdarg.h>
#include <time.h>
#include "../tree/heap.h"

#define DEFAULT_MSG_LEN 10000000
//#include <stdarg.h>

// #define PARALLEL_FIX_BORDERS
// #define WRITE_SHUFFLED_DATA
@@ -60,19 +60,16 @@


#define CHECK_ALLOCATION_NO_CTX(x) if(!x){printf("[!!!] Failed allocation: %s at line %d \n", __FILE__, __LINE__ ); exit(1);}
#ifndef MY_MALLOC
#define MY_MALLOC(n) ({void* p = calloc(n,1); CHECK_ALLOCATION_NO_CTX(p); p; })
#endif

#define DB_PRINT(...) printf(__VA_ARGS__); fflush(stdout)
#ifdef NDEBUG
	#undef DB_PRINT(...)
	#undef DB_PRINT
	#define DB_PRINT(...)
#endif

#define MPI_DB_PRINT(...) mpi_printf(ctx,__VA_ARGS__)
#ifdef NDEBUG
	#undef MPI_DB_PRINT(...)
	#undef MPI_DB_PRINT
	#define MPI_DB_PRINT(...)
#endif

@@ -81,8 +78,8 @@
#ifdef NDEBUG
    #define TIME_DEF 
    #define TIME_START 
    #define TIME_STOP 
    #define LOG_WRITE
    #define TIME_STOP (0.0)
    #define LOG_WRITE(sec_name,time) {}
#else 
    #define TIME_DEF struct timespec __start, __end;
    #define TIME_START { \
@@ -113,22 +110,8 @@
    
#endif

/*
 * from Spriengel code Gadget4
 */

#if defined(__STDC_VERSION__) && __STDC_VERSION__ >= 202000L
/* C2x does not require the second parameter for va_start. */
#define va_start(ap, ...) __builtin_va_start(ap, 0)
#else
/* Versions before C2x do require the second parameter. */
#define va_start(ap, param) __builtin_va_start(ap, param)
#endif
#define va_end(ap)          __builtin_va_end(ap)
#define va_arg(ap, type)    __builtin_va_arg(ap, type)

#if defined(NDEBUG)
    FILE* __log_file;
    static FILE* __log_file;
    #define LOG_START __log_file = fopen("","w"); 
    #define LOG
    #define LOG_END
+73 −60
Original line number Diff line number Diff line
@@ -7,6 +7,7 @@
#include <string.h>
#include <unistd.h>
#include <getopt.h>
#include <errno.h>

#ifdef THREAD_FUNNELED
    #define THREAD_LEVEL MPI_THREAD_FUNNELED
@@ -21,18 +22,18 @@ struct option long_options[] =
    {"in-dtype", required_argument, 0, 't'},
    {"in-dims" , required_argument, 0, 'd'},
    {"out-directory", required_argument, 0, 'o'},
    {"kngbh", optional_argument, 0, 'k'},
    {"z", optional_argument, 0, 'z'},
    {"knn-strategy", optional_argument, 0, 's'},
    {"kngbh", required_argument, 0, 'k'},
    {"z", required_argument, 0, 'z'},
    {"knn-strategy", required_argument, 0, 's'},
    {"parallel-out", no_argument, 0, 'p'},
    {"help", optional_argument, 0, 'h'},
    {"help", no_argument, 0, 'h'},
    {0, 0, 0, 0}
};

const char* help = "Distributed Advanced Density Peak\n"\
                   "    -h --help show this message\n"\
                   "    -i --in-data        (required) path of the input file\n"\
                   "    -t --in-type        (required) datatype of the input file, allowed choices `f32`, `f64`\n"\
                   "    -t --in-dtype       (required) datatype of the input file, allowed choices `f32`, `f64`\n"\
                   "    -d --in-dims        (required) number of dimensions of the data file, dadp expects something\n"\
                   "                                   of the form N x d where N is inferred from the lenght of the\n"\
                   "                                   data file\n"\
@@ -55,7 +56,7 @@ void parse_args(global_context_t* ctx, int argc, char** argv)
        switch(opt)
        {
            case 'i':
                strncpy(ctx -> input_data_file, optarg, DEFAULT_STR_LEN);
                snprintf(ctx -> input_data_file, DEFAULT_STR_LEN, "%s", optarg);
                input_file_set = 1;
                break;
            case 't':
@@ -70,21 +71,38 @@ void parse_args(global_context_t* ctx, int argc, char** argv)
                input_type_set = 1;
                break;
            case 'd':
                ctx -> dims = atoi(optarg);
                if(ctx -> dims < 0)
                {
                    fprintf(stderr, "Invaild number of dimensions\n");
                    char* end;
                    errno = 0;
                    long int dd = strtol(optarg, &end, 10);
                    if(end == optarg || *end != '\0' || errno == ERANGE || dd < 1)
                    {
                        fprintf(stderr, "Invaild number of dimensions (must be >= 1)\n");
                        MPI_Finalize();
                        exit(1);
                    }
                    ctx -> dims = (uint32_t)dd;
                }
                break;
            case 'o':
                strncpy(ctx -> output_data_file, optarg, DEFAULT_STR_LEN);
                snprintf(ctx -> output_data_file, DEFAULT_STR_LEN, "%s", optarg);
                break;
            case 'k':
                ctx -> k = atoi(optarg);
                {
                    char* end;
                    errno = 0;
                    long int kk = strtol(optarg, &end, 10);
                    if(end == optarg || *end != '\0' || errno == ERANGE || kk < 2)
                    {
                        fprintf(stderr, "Invaild number of neighbors (must be >= 2)\n");
                        MPI_Finalize();
                        exit(1);
                    }
                    ctx -> k = (idx_t)kk;
                }
                break;
            case 'z':
                ctx -> z = atof(optarg);
                ctx -> z = strtod(optarg, NULL);
                break;
            case 'h':
                mpi_printf(ctx, "%s\n", help);
@@ -161,7 +179,7 @@ void print_hello(global_context_t* ctx)
}


void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx);
void mpiio_master_read_and_scatter(global_context_t *ctx);

int main(int argc, char** argv) {
    #if defined (_OPENMP)
@@ -202,33 +220,12 @@ int main(int argc, char** argv) {
	ctx.mpi_communicator = MPI_COMM_WORLD;
	get_context(&ctx);


	/*
	 * Mock reading some files, one for each processor
	 */

	int d = 5;
    int k = 300;
	
	float_t* data;
    
    //parse command line
    
    parse_args(&ctx, argc, argv);
    print_hello(&ctx);
	/*
	 * Generate a random matrix of lenght of some kind
	 */


	if(ctx.mpi_rank == 0)
	{
		mpiio_master_read_and_scatter(5, 1000000, &ctx);	
	}
	else
	{
		mpiio_master_read_and_scatter(0, 0, &ctx);	
	}
	mpiio_master_read_and_scatter(&ctx);

	//free(data);
	free_context(&ctx);
@@ -237,7 +234,7 @@ int main(int argc, char** argv) {



void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx) 
void mpiio_master_read_and_scatter(global_context_t *ctx) 
{
    /* TODO
     *
@@ -304,6 +301,25 @@ void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
                                            ctx->input_data_in_float32,
                                            ctx->mpi_communicator, ctx->mpi_rank, ctx->world_size);

    /* input validation that depends on the data (see BUG_REPORT 2.3) */
    {
        size_t min_rank_local_n;
        size_t local_n = (size_t)ctx->local_n_points;
        MPI_Allreduce(&local_n, &min_rank_local_n, 1, MPI_UINT64_T, MPI_MIN, ctx->mpi_communicator);
        if((size_t)ctx->n_points < (size_t)ctx->world_size)
        {
            fprintf(stderr, "n_points (%zu) < world_size (%d): some rank would have no points\n",
                    (size_t)ctx->n_points, ctx->world_size);
            MPI_Abort(ctx->mpi_communicator, 1);
        }
        if((size_t)ctx->k > min_rank_local_n)
        {
            fprintf(stderr, "k (%zu) > min rank local points (%zu): kNN heap would never fill\n",
                    (size_t)ctx->k, min_rank_local_n);
            MPI_Abort(ctx->mpi_communicator, 1);
        }
    }
    
    elapsed_time = TIME_STOP;
    LOG_WRITE("Importing file with MPI-IO", elapsed_time);

@@ -337,19 +353,21 @@ void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)

    TIME_START;

    void* opaque_local_tree;
    union {
        kdtree_t kd;
        vptree_t vp;
    } local_tree;

    switch(ctx -> local_tree_type)
    {
        case KD:
            {
                kdtree_t local_tree;
                kdtree_initialize( &local_tree, ctx -> local_data, ctx -> local_n_points, (unsigned int)ctx -> dims);
                kdtree_initialize( &local_tree.kd, ctx -> local_data, ctx -> local_n_points, (unsigned int)ctx -> dims);


                build_tree_kdtree_parallel(&local_tree);
                // build_tree_kdtree_sample(&local_tree);
                ctx -> local_data = local_tree.data;
                build_tree_kdtree_parallel(&local_tree.kd);
                // build_tree_kdtree_sample(&local_tree.kd);
                ctx -> local_data = local_tree.kd.data;
                
                elapsed_time = TIME_STOP;
                LOG_WRITE("Local kdtrees init and build", elapsed_time);
@@ -358,23 +376,21 @@ void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
                MPI_DB_PRINT("----- Performing ngbh search -----\n");
                MPI_Barrier(ctx -> mpi_communicator);

                mpi_all_knn_search_kdtree(ctx, dp_info, &tree, &local_tree, ctx -> k);
                mpi_all_knn_search_kdtree(ctx, dp_info, &tree, &local_tree.kd, ctx -> k);

                MPI_Barrier(ctx -> mpi_communicator);
                elapsed_time = TIME_STOP;
                LOG_WRITE("Total time for all knn search w. local kdtrees", elapsed_time);
                opaque_local_tree = &local_tree;
            }
            break;
        case VP:
            {
                vptree_t local_tree;
                vptree_initialize( &local_tree, ctx -> local_data, ctx -> local_n_points, (unsigned int)ctx -> dims);
                vptree_initialize( &local_tree.vp, ctx -> local_data, ctx -> local_n_points, (unsigned int)ctx -> dims);


                build_tree_vptree(&local_tree);
                // build_tree_kdtree_sample(&local_tree);
                ctx -> local_data = local_tree.data;
                build_tree_vptree(&local_tree.vp);
                // build_tree_kdtree_sample(&local_tree.vp);
                ctx -> local_data = local_tree.vp.data;

                elapsed_time = TIME_STOP;
                LOG_WRITE("Vptrees init and build", elapsed_time);
@@ -383,12 +399,11 @@ void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
                MPI_DB_PRINT("----- Performing ngbh search -----\n");
                MPI_Barrier(ctx -> mpi_communicator);

                mpi_all_knn_search_vptree(ctx, dp_info, &tree, &local_tree, ctx -> k);
                mpi_all_knn_search_vptree(ctx, dp_info, &tree, &local_tree.vp, ctx -> k);

                MPI_Barrier(ctx -> mpi_communicator);
                elapsed_time = TIME_STOP;
                LOG_WRITE("Total time for all knn search w. vptrees", elapsed_time);
                opaque_local_tree = &local_tree;
            }
            break;
    }
@@ -440,27 +455,25 @@ void mpiio_master_read_and_scatter(int dims, size_t n, global_context_t *ctx)
    {
        case KD:
            {
                kdtree_t local_tree = *((kdtree_t*)opaque_local_tree);
                for(int i = 0; i < ctx -> local_n_points; ++i)
                {
                    cl[i] = ctx -> local_datapoints[i].cluster_idx;
                    idx_t idx = local_tree.__points[i].array_idx;
                    memcpy(data_to_write + idx*ctx -> dims, local_tree.__points[i].data, ctx -> dims*sizeof(float_t));
                    idx_t idx = local_tree.kd.__points[i].array_idx;
                    memcpy(data_to_write + idx*ctx -> dims, local_tree.kd.__points[i].data, ctx -> dims*sizeof(float_t));
                }
                kdtree_free(&local_tree);
                kdtree_free(&local_tree.kd);
            }
            break;

        case VP:
            {
                vptree_t local_tree = *((vptree_t*)opaque_local_tree);
                for(int i = 0; i < ctx -> local_n_points; ++i)
                {
                    cl[i] = ctx -> local_datapoints[i].cluster_idx;
                    idx_t idx = local_tree.__points[i].array_idx;
                    memcpy(data_to_write + idx*ctx -> dims, local_tree.__points[i].data, ctx -> dims*sizeof(float_t));
                    idx_t idx = local_tree.vp.__points[i].array_idx;
                    memcpy(data_to_write + idx*ctx -> dims, local_tree.vp.__points[i].data, ctx -> dims*sizeof(float_t));
                }
                vptree_free(&local_tree);
                vptree_free(&local_tree.vp);
            }
            break;

+19 −1
Original line number Diff line number Diff line
@@ -30,7 +30,25 @@

#define ALIGNMENT 64
#define CHECK_ALLOCATION_NO_CTX(x) if(!x){printf("[!!!] Failed allocation: %s at line %d \n", __FILE__, __LINE__ ); exit(1);}
#define MY_MALLOC(n) ({void* p = aligned_alloc(ALIGNMENT,n); CHECK_ALLOCATION_NO_CTX(p); memset(p, 0, n); p; })
#define MY_MALLOC(n) ({size_t __n = (n + ALIGNMENT - 1) & ~((size_t)ALIGNMENT - 1); void* p = aligned_alloc(ALIGNMENT,__n); CHECK_ALLOCATION_NO_CTX(p); memset(p, 0, __n); p; })

#ifdef _OPENMP
    #include <omp.h>
    static unsigned int __dadp_rng_seed;
    #pragma omp threadprivate(__dadp_rng_seed)
    static inline unsigned int dadp_rand(void)
    {
        if (__dadp_rng_seed == 0)
            __dadp_rng_seed = 42u + (unsigned)omp_get_thread_num();
        return rand_r(&__dadp_rng_seed);
    }
#else
    static unsigned int __dadp_rng_seed = 42u;
    static inline unsigned int dadp_rand(void)
    {
        return rand_r(&__dadp_rng_seed);
    }
#endif


typedef struct {
Loading