Skip to content

Commit

Permalink
localize bucket collection before transfering to main thread
Browse files Browse the repository at this point in the history
  • Loading branch information
@aaronsms
aaronsms committed Dec 5, 2022
1 parent e8f26dd commit 87247c4
Show file tree
Hide file tree
Showing 3 changed files with 136 additions and 77 deletions.
11 changes: 11 additions & 0 deletions raster/r.univar/globals.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,17 @@ typedef struct
int first;
} univar_stat;

/* buffer for each thread to support parallelization */
typedef struct
{
size_t n;
size_t n_alloc;
void *nextp;
DCELL *dcells;
FCELL *fcells;
CELL *cells;
} bucket;

typedef struct
{
CELL min, max, n_zones;
Expand Down
159 changes: 113 additions & 46 deletions raster/r.univar/r.univar_main.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -343,7 +343,7 @@ process_raster(univar_stat * stats, int *fd, int *fdz,

#pragma omp parallel if(nprocs > 1)
{
int i;
int z;
int t_id = 0;

#if defined(_OPENMP)
Expand All @@ -359,9 +359,17 @@ process_raster(univar_stat * stats, int *fd, int *fdz,
double *min = G_malloc(n_alloc * sizeof(double));
double *max = G_malloc(n_alloc * sizeof(double));

for (i = 0; i < n_alloc; i++) {
max[i] = DBL_MIN;
min[i] = DBL_MAX;
bucket *buckets = G_malloc(n_alloc * sizeof(bucket));

for (z = 0; z < n_alloc; z++) {
max[z] = DBL_MIN;
min[z] = DBL_MAX;
buckets[z].n = 0;
buckets[z].n_alloc = 0;
buckets[z].nextp = NULL;
buckets[z].cells = NULL;
buckets[z].fcells = NULL;
buckets[z].dcells = NULL;
}

#pragma omp for schedule(static)
Expand Down Expand Up @@ -404,45 +412,40 @@ process_raster(univar_stat * stats, int *fd, int *fdz,
}

if (param.extended->answer) {
bucket *t_bkt = &buckets[zone];

/* check allocated memory */
/* parallelization is disabled, local variable reflects global state */
if (stats[zone].n + n[zone] >= stats[zone].n_alloc) {
stats[zone].n_alloc += 1000;
if (t_bkt->n >= t_bkt->n_alloc) {
t_bkt->n_alloc += 1000;
size_t msize;

switch (map_type) {
case DCELL_TYPE:
msize = stats[zone].n_alloc * sizeof(DCELL);
stats[zone].dcell_array =
(DCELL *) G_realloc((void *)stats[zone].
dcell_array, msize);
stats[zone].nextp =
(void *)&(stats[zone].dcell_array[stats[zone].n + n[zone]]);
msize = t_bkt->n_alloc * sizeof(DCELL);
t_bkt->dcells =
(DCELL *) G_realloc((void *)t_bkt->dcells,
msize);
t_bkt->nextp = (void *)&(t_bkt->dcells[t_bkt->n]);
break;
case FCELL_TYPE:
msize = stats[zone].n_alloc * sizeof(FCELL);
stats[zone].fcell_array =
(FCELL *) G_realloc((void *)stats[zone].
fcell_array, msize);
stats[zone].nextp =
(void *)&(stats[zone].fcell_array[stats[zone].n + n[zone]]);
msize = t_bkt->n_alloc * sizeof(FCELL);
t_bkt->fcells =
(FCELL *) G_realloc((void *)t_bkt->fcells,
msize);
t_bkt->nextp = (void *)&(t_bkt->fcells[t_bkt->n]);
break;
case CELL_TYPE:
msize = stats[zone].n_alloc * sizeof(CELL);
stats[zone].cell_array =
(CELL *) G_realloc((void *)stats[zone].
cell_array, msize);
stats[zone].nextp =
(void *)&(stats[zone].cell_array[stats[zone].n + n[zone]]);
break;
default:
msize = t_bkt->n_alloc * sizeof(CELL);
t_bkt->cells =
(CELL *) G_realloc((void *)t_bkt->cells,
msize);
t_bkt->nextp = (void *)&(t_bkt->cells[t_bkt->n]);
break;
}
}
/* put the value into stats->XXXcell_array */
memcpy(stats[zone].nextp, ptr, value_sz);
stats[zone].nextp =
G_incr_void_ptr(stats[zone].nextp, value_sz);
memcpy(t_bkt->nextp, ptr, value_sz);
t_bkt->nextp = G_incr_void_ptr(t_bkt->nextp, value_sz);
}

val = ((map_type == DCELL_TYPE) ? *((DCELL *) ptr)
Expand All @@ -461,8 +464,7 @@ process_raster(univar_stat * stats, int *fd, int *fdz,
ptr = G_incr_void_ptr(ptr, value_sz);
if (n_zones)
zptr++;
n[zone]++;

buckets[zone].n++;
} /* end column loop */
if (!(param.shell_style->answer)) {
#pragma omp atomic update
Expand All @@ -471,31 +473,96 @@ process_raster(univar_stat * stats, int *fd, int *fdz,
}
} /* end row loop */

for (i = 0; i < n_alloc; i++) {
for (z = 0; z < n_alloc; z++) {
if (param.extended->answer) {
#pragma omp critical
{
/*
Transfers for each thread from local bucket to global buffer.
Case 1: first transfer, skip reallocation and point to the buffer.
Case 2: other transfers, reallocate exactly if any non-empty bucket.
*/
bucket *t_bucket = &buckets[z];
univar_stat *g_buffer = &stats[z];
size_t old_size;
size_t add_size;

switch (map_type) {
case DCELL_TYPE:
if (NULL == g_buffer->dcell_array) {
g_buffer->dcell_array = t_bucket->dcells;
}
else if (t_bucket->n != 0) {
old_size = g_buffer->n * sizeof(DCELL);
add_size = t_bucket->n * sizeof(DCELL);

g_buffer->dcell_array = (DCELL *) G_realloc((void
*)
g_buffer->dcell_array, old_size + add_size);
memcpy(&g_buffer->dcell_array[g_buffer->n],
t_bucket->dcells, add_size);
}
break;
case FCELL_TYPE:
if (NULL == g_buffer->fcell_array) {
g_buffer->fcell_array = t_bucket->fcells;
}
else if (t_bucket->n != 0) {
old_size = g_buffer->n * sizeof(FCELL);
add_size = t_bucket->n * sizeof(FCELL);

g_buffer->fcell_array = (FCELL *) G_realloc((void
*)
g_buffer->fcell_array, old_size + add_size);
memcpy(&g_buffer->fcell_array[g_buffer->n],
t_bucket->fcells, add_size);
}
break;
case CELL_TYPE:
if (NULL == g_buffer->cell_array) {
g_buffer->cell_array = t_bucket->cells;
}
else if (t_bucket->n != 0) {
old_size = g_buffer->n * sizeof(CELL);
add_size = t_bucket->n * sizeof(CELL);

g_buffer->cell_array = (CELL *) G_realloc((void *)
g_buffer->cell_array, old_size + add_size);
memcpy(&g_buffer->cell_array[g_buffer->n],
t_bucket->cells, add_size);
}
break;
}

g_buffer->n += t_bucket->n;
}
}
else {
#pragma omp atomic update
stats[i].n += n[i];
stats[z].n += buckets[z].n;
}
#pragma omp atomic update
stats[i].size += size[i];
stats[z].size += size[z];
#pragma omp atomic update
stats[i].sum += sum[i];
stats[z].sum += sum[z];
#pragma omp atomic update
stats[i].sumsq += sumsq[i];
stats[z].sumsq += sumsq[z];
#pragma omp atomic update
stats[i].sum_abs += sum_abs[i];
stats[z].sum_abs += sum_abs[z];

#if defined(_OPENMP)
omp_set_lock(&(minmax[i]));
omp_set_lock(&(minmax[z]));
#endif
if (stats[i].max < max[i] ||
(stats[i].max != stats[i].max && max[i] != DBL_MIN)) {
stats[i].max = max[i];
if (stats[z].max < max[z] ||
(stats[z].max != stats[z].max && max[z] != DBL_MIN)) {
stats[z].max = max[z];
}
if (stats[i].min > min[i] ||
(stats[i].min != stats[i].min && min[i] != DBL_MAX)) {
stats[i].min = min[i];
if (stats[z].min > min[z] ||
(stats[z].min != stats[z].min && min[z] != DBL_MAX)) {
stats[z].min = min[z];
}
#if defined(_OPENMP)
omp_unset_lock(&(minmax[i]));
omp_unset_lock(&(minmax[z]));
#endif
}

Expand Down
43 changes: 12 additions & 31 deletions raster/r.univar/stats.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,27 +44,8 @@ univar_stat *create_univar_stat_struct(int map_type, int n_perc)
stats[i].fcell_array = NULL;
stats[i].cell_array = NULL;
stats[i].map_type = map_type;

stats[i].n_alloc = 0;

stats[i].first = TRUE;

/* allocate memory for extended computation */
/* changed to on-demand block allocation */

/* if (param.extended->answer) {
if (map_type == DCELL_TYPE) {
stats[i].dcell_array = NULL;
}
else if (map_type == FCELL_TYPE) {
stats[i].fcell_array =NULL;
}
else if (map_type == CELL_TYPE) {
stats[i].cell_array = NULL;
}
}
*/

}

return stats;
Expand Down Expand Up @@ -212,8 +193,8 @@ int print_stats(univar_stat * stats)
quartile_25 = (double)stats[z].cell_array[qpos_25];
if (stats[z].n % 2) /* odd */
median =
(double)stats[z].
cell_array[(int)(stats[z].n / 2)];
(double)
stats[z].cell_array[(int)(stats[z].n / 2)];
else /* even */
median =
(double)(stats[z].cell_array[stats[z].n / 2 - 1] +
Expand All @@ -232,12 +213,12 @@ int print_stats(univar_stat * stats)
quartile_25 = (double)stats[z].fcell_array[qpos_25];
if (stats[z].n % 2) /* odd */
median =
(double)stats[z].
fcell_array[(int)(stats[z].n / 2)];
(double)
stats[z].fcell_array[(int)(stats[z].n / 2)];
else /* even */
median =
(double)(stats[z].
fcell_array[stats[z].n / 2 - 1] +
(double)(stats[z].fcell_array[stats[z].n / 2 - 1]
+
stats[z].fcell_array[stats[z].n / 2]) /
2.0;
quartile_75 = (double)stats[z].fcell_array[qpos_75];
Expand Down Expand Up @@ -478,8 +459,8 @@ int print_stats_table(univar_stat * stats)
quartile_25 = (double)stats[z].cell_array[qpos_25];
if (stats[z].n % 2) /* odd */
median =
(double)stats[z].
cell_array[(int)(stats[z].n / 2)];
(double)
stats[z].cell_array[(int)(stats[z].n / 2)];
else /* even */
median =
(double)(stats[z].cell_array[stats[z].n / 2 - 1] +
Expand All @@ -498,12 +479,12 @@ int print_stats_table(univar_stat * stats)
quartile_25 = (double)stats[z].fcell_array[qpos_25];
if (stats[z].n % 2) /* odd */
median =
(double)stats[z].
fcell_array[(int)(stats[z].n / 2)];
(double)
stats[z].fcell_array[(int)(stats[z].n / 2)];
else /* even */
median =
(double)(stats[z].
fcell_array[stats[z].n / 2 - 1] +
(double)(stats[z].fcell_array[stats[z].n / 2 - 1]
+
stats[z].fcell_array[stats[z].n / 2]) /
2.0;
quartile_75 = (double)stats[z].fcell_array[qpos_75];
Expand Down

0 comments on commit 87247c4

Please sign in to comment.