stop program if invalid array size occurs
[pub/astra/parallel.git] / prefix / hillis_sum.c
1 /*
2  * O(nlog n) work algorithm (Hillis-Steele)
3  */
4 #include "hillis_sum.h"
5
6 /*
7  * Hillis/Steele, prefix sum version
8  */
9 unsigned long is_power_of_two(unsigned long var) {
10         if (var == 1) {
11                 return 1;
12         }
13         while (var > 2) {
14                 if ((var & 0x1) == 0x1) {
15                         return 0;
16                 }
17                 var >>= 1;
18         }
19         return 1;
20 }
21
22 void algorithm (numtype x[], unsigned long size, unsigned int ops[]) {
23         unsigned long k;
24         unsigned long i;
25
26         if (is_power_of_two(size) != 1) {
27                 fprintf(stderr, "Hillis-Sum algorithm only works on arrays on size 2^n\n");
28                 exit(1);
29         }
30
31         for(k=2; k <= size; k <<=1){
32                 #pragma omp parallel for shared(x, size, ops, k) private(i)
33                 for(i = (k-1); i < size; i+=k){
34 #ifdef DEBUG
35                         printf ("x[%2li] = x[%2li] + x[%2li]; // {i:%li, k:%li}\n", i, i-k, i, i, k);
36 #endif
37                         x[i] = x[i-(k/2)] + x[i];
38                 }
39         }
40 }
41