/*Scalable Universal Matrix Multiplication Algorithm) works by dividing matrices into blocks and processing them in phases.*/
|
|
/*Matrix size n is divisible by sqrt(p) where p = number of threads.*/
|
|
/*Flat (1D) memory is used for all matrices.*/
|
|
/*Blocked computation is handled with explicit loop tiling. -1990*/
|
|
/*Partitions the matrices A and B across a 2D processor grid
|
|
|
|
Performs a series of broadcasts of matrix blocks along rows and columns
|
|
|
|
Computes partial results at each processor
|
|
|
|
Accumulates results to form the final product matrix C = A × B*/
|
|
#include <iostream>
|
|
#include <cstdlib>
|
|
#include <cmath>
|
|
#include <random>
|
|
#include <boost/thread.hpp>
|
|
#include <boost/bind/bind.hpp>
|
|
#include <boost/chrono.hpp>
|
|
|
|
using namespace std;
|
|
using namespace boost::placeholders;
|
|
|
|
const int MAX_THREADS = 100;
|
|
|
|
// Allocate 1D array matrix
|
|
double* create_matrix(int n) {
|
|
return new double[n * n];
|
|
}
|
|
|
|
void cleanup_matrix(double* matrix) {
|
|
delete[] matrix;
|
|
}
|
|
|
|
void init_matrix(double* matrix, int n, bool zero = false) {
|
|
random_device rd;
|
|
mt19937 gen(rd());
|
|
uniform_real_distribution<double> dis(0.0, 1.0);
|
|
|
|
for (int i = 0; i < n * n; ++i)
|
|
matrix[i] = zero ? 0.0 : dis(gen);
|
|
}
|
|
|
|
// SUMMA worker thread: computes block (i_block) of C
|
|
void summa_block(int block_row_start, int block_row_end, int n, int block_size,
|
|
double* a, double* b, double* c) {
|
|
for (int i = block_row_start; i < block_row_end; ++i) {
|
|
for (int j = 0; j < n; ++j) {
|
|
double sum = 0.0;
|
|
for (int k = 0; k < n; ++k) {
|
|
sum += a[i * n + k] * b[k * n + j];
|
|
}
|
|
c[i * n + j] = sum;
|
|
}
|
|
}
|
|
}
|
|
|
|
int main(int argc, char** argv) {
|
|
if (argc != 2) {
|
|
cerr << "Usage: " << argv[0] << " <matrix_size>" << endl;
|
|
return 1;
|
|
}
|
|
|
|
int n = atoi(argv[1]);
|
|
if (n <= 0) {
|
|
cerr << "Matrix size must be positive" << endl;
|
|
return 1;
|
|
}
|
|
|
|
boost::chrono::system_clock::time_point t1, t2, t3, t4;
|
|
t1 = boost::chrono::system_clock::now();
|
|
|
|
double* a = create_matrix(n);
|
|
double* b = create_matrix(n);
|
|
double* c = create_matrix(n);
|
|
|
|
init_matrix(a, n);
|
|
init_matrix(b, n);
|
|
init_matrix(c, n, true);
|
|
|
|
t2 = boost::chrono::system_clock::now();
|
|
t3 = boost::chrono::system_clock::now();
|
|
|
|
boost::thread_group group;
|
|
|
|
int rows_per_thread = n / MAX_THREADS;
|
|
int remaining = n % MAX_THREADS;
|
|
int current_row = 0;
|
|
|
|
for (int i = 0; i < MAX_THREADS && current_row < n; ++i) {
|
|
int rows_to_process = rows_per_thread + (i < remaining ? 1 : 0);
|
|
int end_row = current_row + rows_to_process;
|
|
|
|
group.create_thread(
|
|
boost::bind(summa_block, current_row, end_row, n, rows_to_process, a, b, c)
|
|
);
|
|
|
|
current_row = end_row;
|
|
}
|
|
|
|
group.join_all();
|
|
|
|
t4 = boost::chrono::system_clock::now();
|
|
|
|
auto setup_time = boost::chrono::duration_cast<boost::chrono::milliseconds>(t2 - t1);
|
|
auto exec_time = boost::chrono::duration_cast<boost::chrono::milliseconds>(t4 - t3);
|
|
|
|
cout << "n,setup(ms),exec(ms)" << endl;
|
|
cout << n << "," << setup_time.count() << "," << exec_time.count() << endl;
|
|
|
|
cleanup_matrix(a);
|
|
cleanup_matrix(b);
|
|
cleanup_matrix(c);
|
|
|
|
return 0;
|
|
}
|
|
|