13 #include "../offload/offload_runtime.h"
26 static inline int imax(
int x,
int y) {
return (x > y ? x : y); }
32 static inline void min_max(
int result[2],
int value) {
33 if (value < result[0]) {
36 if (result[1] < value) {
46 const double filter_eps) {
47 const int nrows = (trans) ? matrix->
ncols : matrix->
nrows;
48 int *nblocks_per_row = calloc(nrows,
sizeof(
int));
49 float *row_max_eps = malloc(nrows *
sizeof(
float));
56 for (
int iblock = 0; iblock < shard->
nblocks; iblock++) {
58 const int row = (trans) ? blk->
col : blk->
row;
60 nblocks_per_row[row]++;
67 for (
int i = 0;
i < nrows;
i++) {
69 ((float)filter_eps) / ((float)
imax(1, nblocks_per_row[
i]));
70 row_max_eps[
i] = f * f;
74 free(nblocks_per_row);
83 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
84 dbm_multiply_gpu_context_t gpu;
95 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
97 matrix_c->
shards, &ctx->gpu);
113 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
114 dbm_multiply_gpu_upload_packs(pack_a, pack_b, &ctx->gpu);
127 const int mnk_range[3][2],
const double alpha,
132 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
136 dbm_multiply_gpu_process_batch(ntasks, batch, mnk_range, alpha, kshard,
151 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
152 dbm_multiply_gpu_download_results(&ctx->gpu);
163 #if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_DBM)
164 dbm_multiply_gpu_stop(&ctx->gpu);
178 const bool retain_sparsity,
179 const float *rows_max_eps, int64_t *flop,
181 const float alpha2 = alpha * alpha;
182 int64_t flop_sum = 0;
186 int shard_row_start[nshard_rows], shard_col_start[nshard_cols];
187 memset(shard_row_start, 0, nshard_rows *
sizeof(
int));
188 memset(shard_col_start, 0, nshard_cols *
sizeof(
int));
190 const int *sum_index_sizes_a =
192 const int *sum_index_sizes_b =
194 const int *free_index_sizes_a =
196 const int *free_index_sizes_b =
199 #pragma omp parallel reduction(+ : flop_sum)
204 for (
int iblock = 1; iblock < pack_a->
nblocks; iblock++) {
206 const int prev_shard_row =
208 if (prev_shard_row != shard_row) {
209 shard_row_start[shard_row] = iblock;
213 for (
int jblock = 1; jblock < pack_b->
nblocks; jblock++) {
215 const int prev_shard_col =
217 if (prev_shard_col != shard_col) {
218 shard_col_start[shard_col] = jblock;
222 #pragma omp for collapse(2) schedule(dynamic)
223 for (
int shard_row = 0; shard_row < nshard_rows; shard_row++) {
224 for (
int shard_col = 0; shard_col < nshard_cols; shard_col++) {
225 const int ishard = shard_row * nshard_cols + shard_col;
228 int mnk_range[][2] = {{INT_MAX, 0}, {INT_MAX, 0}, {INT_MAX, 0}};
233 const int iblock_start = shard_row_start[shard_row];
234 int jblock_start = shard_col_start[shard_col];
235 for (
int iblock = iblock_start; iblock < pack_a->
nblocks; iblock++) {
237 if (blk_a->
free_index % nshard_rows != shard_row) {
240 for (
int jblock = jblock_start; jblock < pack_b->
nblocks; jblock++) {
242 if (blk_b->
free_index % nshard_cols != shard_col) {
255 const float result_norm = alpha2 * blk_a->
norm * blk_b->
norm;
256 if (result_norm < rows_max_eps[blk_a->
free_index]) {
261 const int m = free_index_sizes_a[blk_a->
free_index];
262 const int n = free_index_sizes_b[blk_b->
free_index];
263 const int k = sum_index_sizes_a[blk_a->
sum_index];
266 assert(k == sum_index_sizes_b[blk_b->
sum_index]);
271 if (blk_c == NULL && retain_sparsity) {
273 }
else if (blk_c == NULL) {
282 const int task_flops = 2 * m * n * k;
283 flop_sum += task_flops;
284 if (task_flops == 0) {
304 pack_b, ishard, shard_c, ctx);
305 mnk_range[0][0] = mnk_range[1][0] = mnk_range[2][0] = INT_MAX;
306 mnk_range[0][1] = mnk_range[1][1] = mnk_range[2][1] = 0;
312 ishard, shard_c, ctx);
324 void dbm_multiply(
const bool transa,
const bool transb,
const double alpha,
327 const bool retain_sparsity,
const double filter_eps,
330 assert(omp_get_num_threads() == 1);
334 const int num_sum_index_a = (transa) ? matrix_a->
nrows : matrix_a->
ncols;
335 const int num_sum_index_b = (transb) ? matrix_b->
ncols : matrix_b->
nrows;
336 const int num_free_index_a = (transa) ? matrix_a->
ncols : matrix_a->
nrows;
337 const int num_free_index_b = (transb) ? matrix_b->
nrows : matrix_b->
ncols;
340 assert(num_sum_index_a == num_sum_index_b);
341 assert(num_free_index_a == matrix_c->
nrows);
342 assert(num_free_index_b == matrix_c->
ncols);
362 multiply_packs(transa, transb, alpha, pack_a, pack_b, matrix_a, matrix_b,
363 matrix_c, retain_sparsity, rows_max_eps, flop, ctx);
static const int MAX_BATCH_SIZE
void dbm_library_counter_increment(const int m, const int n, const int k)
Add given block multiplication to stats. This routine is thread-safe.
static int dbm_get_shard_index(const dbm_matrix_t *matrix, const int row, const int col)
Internal routine for getting a block's shard index.
static int dbm_get_num_shards(const dbm_matrix_t *matrix)
Internal routine that returns the number of shards for given matrix.
void dbm_mpi_sum_int64(int64_t *values, const int count, const dbm_mpi_comm_t comm)
Wrapper around MPI_Allreduce for op MPI_SUM and datatype MPI_INT64_T.
void dbm_mpi_sum_int(int *values, const int count, const dbm_mpi_comm_t comm)
Wrapper around MPI_Allreduce for op MPI_SUM and datatype MPI_INT.
static void backend_upload_packs(const dbm_pack_t *pack_a, const dbm_pack_t *pack_b, backend_context_t *ctx)
Private routine for handing newly arrived packs to the backend.
static void min_max(int result[2], int value)
Updates the min/max of a range of values (initially {INT_MAX, 0}).
static float * compute_rows_max_eps(const bool trans, const dbm_matrix_t *matrix, const double filter_eps)
Private routine for computing the max filter threshold for each row.
static void backend_process_batch(const int ntasks, dbm_task_t batch[ntasks], const int mnk_range[3][2], const double alpha, const dbm_pack_t *pack_a, const dbm_pack_t *pack_b, const int kshard, dbm_shard_t *shard_c, backend_context_t *ctx)
Private routine for sending a batch to the multiplication backend.
static void backend_stop(backend_context_t *ctx)
Private routine for shutting down the multiplication backend.
static void backend_download_results(backend_context_t *ctx)
Private routine for downloading results of the multiplication backend.
static backend_context_t * backend_start(const dbm_matrix_t *matrix_c)
Private routine for intializing the multiplication backend.
static int imax(int x, int y)
Returns the larger of two given integer (missing from the C standard).
static void multiply_packs(const bool transa, const bool transb, const double alpha, const dbm_pack_t *pack_a, const dbm_pack_t *pack_b, const dbm_matrix_t *matrix_a, const dbm_matrix_t *matrix_b, dbm_matrix_t *matrix_c, const bool retain_sparsity, const float *rows_max_eps, int64_t *flop, backend_context_t *ctx)
Private routine for multipling two packs.
void dbm_multiply(const bool transa, const bool transb, const double alpha, const dbm_matrix_t *matrix_a, const dbm_matrix_t *matrix_b, const double beta, dbm_matrix_t *matrix_c, const bool retain_sparsity, const double filter_eps, int64_t *flop)
Performs a multiplication of two dbm_matrix_t matrices. See dbm_matrix.h for details.
dbm_comm_iterator_t * dbm_comm_iterator_start(const bool transa, const bool transb, const dbm_matrix_t *matrix_a, const dbm_matrix_t *matrix_b, const dbm_matrix_t *matrix_c)
Internal routine for creating a communication iterator.
void dbm_comm_iterator_stop(dbm_comm_iterator_t *iter)
Internal routine for releasing the given communication iterator.
bool dbm_comm_iterator_next(dbm_comm_iterator_t *iter, dbm_pack_t **pack_a, dbm_pack_t **pack_b)
Internal routine for retriving next pair of packs from given iterator.
void dbm_multiply_cpu_process_batch(const int ntasks, dbm_task_t batch[ntasks], const double alpha, const dbm_pack_t *pack_a, const dbm_pack_t *pack_b, dbm_shard_t *shard_c)
Internal routine for executing the tasks in given batch on the CPU.
dbm_block_t * dbm_shard_promise_new_block(dbm_shard_t *shard, const int row, const int col, const int block_size)
Internal routine for allocating the metadata of a new block.
dbm_block_t * dbm_shard_lookup(const dbm_shard_t *shard, const int row, const int col)
Internal routine for looking up a block from a shard.
static void const int const int i
subroutine, public dbm_scale(matrix, alpha)
Multiplies all entries in the given matrix by the given factor alpha.
subroutine, public dbm_filter(matrix, eps)
Removes all blocks from the given matrix whose block norm is below the given threshold....
subroutine, public dbm_get_stored_coordinates(matrix, row, column, processor)
Returns the MPI rank on which the given block should be stored.
Private struct for storing the context of the multiplication backend.
Internal struct for storing a block's metadata.
Internal struct for storing a communication iterator.
Internal struct for storing a matrix.
dbm_distribution_t * dist
Internal struct for storing a dbm_block_t plus its norm.
Internal struct for storing a pack - essentially a shard for MPI.
dbm_pack_block_t * blocks
Internal struct for storing a matrix shard.
Internal struct for storing a task, ie. a single block multiplication.