#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdint.h>
#include <omp.h>


size_t get_subtree_size(const size_t array_size) {                                  // calculates subtree size
    size_t subtree_size = 1;
    while (subtree_size < array_size) {
        subtree_size <<= 1u;
    }
    return subtree_size;
}


uint64_t * initialize_array(const size_t array_size) {                              // initialize array
    uint64_t *a = (uint64_t *) calloc(array_size, sizeof(uint64_t));

    for (size_t i = 0; i < array_size; ++i) {
        a[i] = i;
    }

    return a;
}


void build_range_tree(
    const uint64_t *a, const size_t array_size, uint64_t *tree,
    const size_t subtree_size, const size_t threads_count) {

    #pragma omp parallel num_threads(threads_count)
    {
        #pragma omp for schedule(static)
        for (size_t i = 0; i < array_size; ++i) {
            tree[i + subtree_size] = a[i];
        }
    }

    size_t current_layer_size = (subtree_size >> 1u);
    while (current_layer_size > 0) {
        #pragma opm parallel num_threads(threads_count)
        {
            #pragma omp for schedule(static)
            for (size_t i = current_layer_size; i < (current_layer_size << 1u); ++i) {
                tree[i] = tree[i << 1u] + tree[(i << 1u) + 1];
            }
        }
        current_layer_size >>= 1u;
    }
}


uint64_t find_sum_from_subtrees(                                               // right boundaries are not included
    const uint64_t* tree, const size_t left, const size_t right,
    const size_t current_node_index, const size_t current_left,
    const size_t current_right) {

    if (left == right || right < current_left || current_right < left) {       // target range not in current observed range
        return 0;
    }
    if (current_left <= left && right <= current_right) {                      // target range is completeley in observed range
        return tree[current_node_index];
    }

    const size_t middle = (right + left) >> 1u;                                // getting sums from subtrees
    return
        find_sum_from_subtrees(tree, left, right, current_node_index << 1u, current_left, middle) +
        find_sum_from_subtrees(tree, left, right, (current_node_index << 1u) + 1, middle, current_right);
}


uint64_t sum_query(                                                            // query wrapper function (for convenient usage)
    const uint64_t *tree, const size_t subtree_size,
    const size_t left, const size_t right) {                                   // right boundary is included

    return find_sum_from_subtrees(tree, left, right + 1, 1, 0, subtree_size);
}


void perform_queries(                                                          // function for performing queries
    const uint64_t *tree, const size_t queries_count,
    uint64_t *results, const size_t subtree_size, const size_t threads_count) {

    #pragma omp parallel num_threads(threads_count)
    {
        #pragma omp for schedule(static)
        for (size_t i = 0; i < queries_count; ++i) {
            size_t left_bound = (subtree_size >> 2u);
            size_t right_bound = left_bound + (subtree_size >> 1u);
            results[i] = sum_query(tree, subtree_size, left_bound, right_bound);
        }
    }
}


int main(const int argc, const char** argv) {
    if (argc != 5) {                                // check valid launch parameters
        printf(
            "Invalid arguments!\n"
            "argv[1] = threads count\n"
            "argv[2] = input array size\n"
            "argv[3] = queries count\n"
            "argv[4] = launches' count\n"
        );
        return 1;
    }

    int32_t threads_count;
    size_t array_size;
    size_t queries_count;
    size_t launches_count;
    size_t subtree_size;

    sscanf(argv[1], "%d", &threads_count);         // reading program parameters
    if (threads_count <= 0) {
        printf("threads count must be > 0");
    return 2;
    }
    sscanf(argv[2], "%llu", &array_size);
    sscanf(argv[3], "%llu", &queries_count);
    sscanf(argv[4], "%llu", &launches_count);
    subtree_size = get_subtree_size(array_size);

    uint64_t *a = initialize_array(array_size);                                 // array initialization
    uint64_t *results = (uint64_t *) calloc(queries_count, sizeof(uint64_t));   // getting memroy for calculation results holding
    uint64_t *tree = (uint64_t *) calloc(subtree_size << 1u, sizeof(uint64_t)); // getting memory for tree holding

    double summary_building_time = 0;
    double summary_queries_time = 0;

    for (size_t i = 0; i < launches_count; ++i) {                                   // main launch loop for time measurement
        double building_start = omp_get_wtime();                                    // building time measurement
        build_range_tree(a, array_size, tree, subtree_size, (size_t) threads_count);
    summary_building_time += (omp_get_wtime() - building_start);

    double queries_start = omp_get_wtime();                                             // queries execution time measurement
        perform_queries(tree, queries_count, results, subtree_size, (size_t) threads_count);
    summary_queries_time += (omp_get_wtime() - queries_start);
    }

    free(a);
    free(results);
    free(tree);

    printf("%lf %lf\n", summary_building_time / launches_count, summary_queries_time / launches_count);
    return 0;
}