diff --git a/src/serial _real.c b/src/serial _real.c
new file mode 100644
index 0000000000000000000000000000000000000000..ab02e0f5a4eb458acf82d82df7410930d60ba64f
--- /dev/null
+++ b/src/serial _real.c
@@ -0,0 +1,277 @@
+// serial.c
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <time.h>
+
+#define NMAX 100
+#define DATAMAX 1000
+#define DATAMIN -1000
+
+/*
+ * Struct Matrix
+ *
+ * Matrix representation consists of matrix data
+ * and effective dimensions
+ * */
+typedef struct Matrix {
+ int mat[NMAX][NMAX]; // Matrix cells
+ int row_eff; // Matrix effective row
+ int col_eff; // Matrix effective column
+} Matrix;
+
+
+/*
+ * Procedure init_matrix
+ *
+ * Initializing newly allocated matrix
+ * Setting all data to 0 and effective dimensions according
+ * to nrow and ncol
+ * */
+void init_matrix(Matrix *m, int nrow, int ncol) {
+ m->row_eff = nrow;
+ m->col_eff = ncol;
+
+ for (int i = 0; i < m->row_eff; i++) {
+ for (int j = 0; j < m->col_eff; j++) {
+ m->mat[i][j] = 0;
+ }
+ }
+}
+
+
+/*
+ * Function input_matrix
+ *
+ * Returns a matrix with values from stdin input
+ * */
+Matrix input_matrix(int nrow, int ncol) {
+ Matrix input;
+ init_matrix(&input, nrow, ncol);
+
+ for (int i = 0; i < nrow; i++) {
+ for (int j = 0; j < ncol; j++) {
+ scanf("%d", &input.mat[i][j]);
+ }
+ }
+
+ return input;
+}
+
+
+/*
+ * Procedure print_matrix
+ *
+ * Print matrix data
+ * */
+void print_matrix(Matrix *m) {
+ for (int i = 0; i < m->row_eff; i++) {
+ for (int j = 0; j < m->col_eff; j++) {
+ printf("%d ", m->mat[i][j]);
+ }
+ printf("\n");
+ }
+}
+
+
+/*
+ * Function get_matrix_datarange
+ *
+ * Returns the range between maximum and minimum
+ * element of a matrix
+ * */
+int get_matrix_datarange(Matrix *m) {
+ int max = DATAMIN;
+ int min = DATAMAX;
+ for (int i = 0; i < m->row_eff; i++) {
+ for (int j = 0; j < m->col_eff; j++) {
+ int el = m->mat[i][j];
+ if (el > max) max = el;
+ if (el < min) min = el;
+ }
+ }
+
+ return max - min;
+}
+
+
+/*
+ * Function supression_op
+ *
+ * Returns the sum of intermediate value of special multiplication
+ * operation where kernel[0][0] corresponds to target[row][col]
+ * */
+int supression_op(Matrix *kernel, Matrix *target, int row, int col) {
+ int intermediate_sum = 0;
+ for (int i = 0; i < kernel->row_eff; i++) {
+ for (int j = 0; j < kernel->col_eff; j++) {
+ intermediate_sum += kernel->mat[i][j] * target->mat[row + i][col + j];
+ }
+ }
+
+ return intermediate_sum;
+}
+
+
+/*
+ * Function convolution
+ *
+ * Return the output matrix of convolution operation
+ * between kernel and target
+ * */
+Matrix convolution(Matrix *kernel, Matrix *target) {
+ Matrix out;
+ int out_row_eff = target->row_eff - kernel->row_eff + 1;
+ int out_col_eff = target->col_eff - kernel->col_eff + 1;
+
+ init_matrix(&out, out_row_eff, out_col_eff);
+
+ for (int i = 0; i < out.row_eff; i++) {
+ for (int j = 0; j < out.col_eff; j++) {
+ out.mat[i][j] = supression_op(kernel, target, i, j);
+ }
+ }
+
+ return out;
+}
+
+
+/*
+ * Procedure merge_array
+ *
+ * Merges two subarrays of n with n[left..mid] and n[mid+1..right]
+ * to n itself, with n now ordered ascendingly
+ * */
+void merge_array(int *n, int left, int mid, int right) {
+ int n_left = mid - left + 1;
+ int n_right = right - mid;
+ int iter_left = 0, iter_right = 0, iter_merged = left;
+ int arr_left[n_left], arr_right[n_right];
+
+ for (int i = 0; i < n_left; i++) {
+ arr_left[i] = n[i + left];
+ }
+
+ for (int i = 0; i < n_right; i++) {
+ arr_right[i] = n[i + mid + 1];
+ }
+
+ while (iter_left < n_left && iter_right < n_right) {
+ if (arr_left[iter_left] <= arr_right[iter_right]) {
+ n[iter_merged] = arr_left[iter_left++];
+ } else {
+ n[iter_merged] = arr_right[iter_right++];
+ }
+ iter_merged++;
+ }
+
+ while (iter_left < n_left) {
+ n[iter_merged++] = arr_left[iter_left++];
+ }
+ while (iter_right < n_right) {
+ n[iter_merged++] = arr_right[iter_right++];
+ }
+}
+
+
+/*
+ * Procedure merge_sort
+ *
+ * Sorts array n with merge sort algorithm
+ * */
+void merge_sort(int *n, int left, int right) {
+ if (left < right) {
+ int mid = left + (right - left) / 2;
+
+ merge_sort(n, left, mid);
+ merge_sort(n, mid + 1, right);
+
+ merge_array(n, left, mid, right);
+ }
+}
+
+
+/*
+ * Procedure print_array
+ *
+ * Prints all elements of array n of size to stdout
+ * */
+void print_array(int *n, int size) {
+ for (int i = 0; i < size; i++ ) printf("%d ", n[i]);
+ printf("\n");
+}
+
+
+/*
+ * Function get_median
+ *
+ * Returns median of array n of length
+ * */
+int get_median(int *n, int length) {
+ int mid = length / 2;
+ if (length & 1) return n[mid];
+
+ return (n[mid - 1] + n[mid]) / 2;
+}
+
+
+/*
+ * Function get_floored_mean
+ *
+ * Returns floored mean from an array of integers
+ * */
+long get_floored_mean(int *n, int length) {
+ long sum = 0;
+ for (int i = 0; i < length; i++) {
+ sum += n[i];
+ }
+
+ return sum / length;
+}
+
+
+
+// main() driver
+int main() {
+ // Time.
+ clock_t t;
+ t = clock();
+
+ int kernel_row, kernel_col, target_row, target_col, num_targets;
+
+ // reads kernel's row and column and initalize kernel matrix from input
+ scanf("%d %d", &kernel_row, &kernel_col);
+ Matrix kernel = input_matrix(kernel_row, kernel_col);
+
+ // reads number of target matrices and their dimensions.
+ // initialize array of matrices and array of data ranges (int)
+ scanf("%d %d %d", &num_targets, &target_row, &target_col);
+ Matrix* arr_mat = (Matrix*)malloc(num_targets * sizeof(Matrix));
+ int arr_range[num_targets];
+
+ // read each target matrix, compute their convolution matrices, and compute their data ranges
+ for (int i = 0; i < num_targets; i++) {
+ arr_mat[i] = input_matrix(target_row, target_col);
+ arr_mat[i] = convolution(&kernel, &arr_mat[i]);
+ arr_range[i] = get_matrix_datarange(&arr_mat[i]);
+ }
+
+ // sort the data range array
+ merge_sort(arr_range, 0, num_targets - 1);
+
+ int median = get_median(arr_range, num_targets);
+ int floored_mean = get_floored_mean(arr_range, num_targets);
+
+ // print the min, max, median, and floored mean of data range array
+ printf("%d\n%d\n%d\n%d\n",
+ arr_range[0],
+ arr_range[num_targets - 1],
+ median,
+ floored_mean);
+
+ // Print execution time in seconds.
+ t = clock() - t;
+ printf("\nRuntime: %f s\n", ((float)t) / CLOCKS_PER_SEC);
+
+ return 0;
+}
diff --git a/src/serial.c b/src/serial.c
index ab02e0f5a4eb458acf82d82df7410930d60ba64f..f4a69bd39472d4c3c9519557833c2a6e7981e22b 100644
--- a/src/serial.c
+++ b/src/serial.c
@@ -3,6 +3,7 @@
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
+#include <math.h>
#define NMAX 100
#define DATAMAX 1000
@@ -229,6 +230,134 @@ long get_floored_mean(int *n, int length) {
return sum / length;
}
+/**
+ * Function index_to_row_major
+ *
+ * Returns the index of a matrix element in row-major order
+ */
+int index_to_row_major(int row, int col, int row_eff, int col_eff) {
+ return row * col_eff + col;
+}
+
+// __device__ int d_index_to_row_major(int row, int col, int row_eff, int col_eff) {
+// return row * col_eff + col;
+// }
+
+/**
+ * Function row_major_to_index
+ *
+ * Returns the row and column of a matrix element in row-major order
+ */
+void row_major_to_index(int index, int row_eff, int col_eff, int *row, int *col) {
+ *row = index / col_eff;
+ *col = index % col_eff;
+}
+
+// __device__ void d_row_major_to_index(int index, int row_eff, int col_eff, int *row, int *col) {
+// *row = index / col_eff;
+// *col = index % col_eff;
+// }
+
+/**
+ * Function map_matrix
+ *
+ * Returns a row major matrix of the input matrix.
+ **/
+int* map_matrix(int mat[][100], int row, int col) {
+ int* map = (int*) malloc(sizeof(int) * row * col);
+ for (int i = 0; i < row; i++) {
+ for (int j = 0; j < col; j++) {
+ map[index_to_row_major(i, j, row, col)] = mat[i][j];
+ }
+ }
+ return map;
+}
+
+/**
+ * Function reverse_map_matrix
+ *
+ * Returns a matrix of the input row major matrix.
+ */
+int** reverse_map_matrix(int* map, int row, int col) {
+ int** mat = (int**) malloc(sizeof(int*) * row);
+ for (int i = 0; i < row; i++) {
+ mat[i] = (int*) malloc(sizeof(int) * col);
+ for (int j = 0; j < col; j++) {
+ mat[i][j] = map[index_to_row_major(i, j, row, col)];
+ }
+ }
+ return mat;
+}
+
+/**
+ * Function rm_to_matrix_object
+ *
+ * Return Matrix struct of row major matrix
+ */
+Matrix rm_to_matrix_object(int* map, int row, int col) {
+ Matrix mat;
+ init_matrix(&mat, row, col);
+ for (int i = 0; i < row; i++) {
+ for (int j = 0; j < col; j++) {
+ mat.mat[i][j] = map[index_to_row_major(i, j, row, col)];
+ }
+ }
+ return mat;
+}
+
+/**
+ * Function cuda_convolution
+ *
+ * Returns a matrix of the convolution of the input matrix with the kernel
+ */
+void cuda_convolution(int* d_out_mat, int* arr_mat_rm, int* kernel_rm, int row_eff, int col_eff, int kernel_row, int kernel_col) {
+ // Calculate real row and column of input matrix.
+ int row = row_eff + kernel_row - 1;
+ int col = col_eff + kernel_col - 1;
+
+ // For each element in input matrix that is not on the boundary,
+ for (int i = 0; i < row_eff; i++) {
+ for (int j = 0; j < col_eff; j++) {
+ // Convolution of the element with the kernel.
+ // Calculate the sum of the kernel and the input matrix.
+ int intermediate_sum = 0;
+ for (int k = 0; k < kernel_row; k++) {
+ for (int l = 0; l < kernel_col; l++) {
+ int index = index_to_row_major(i + k, j + l, row, col);
+ int kernel_index = index_to_row_major(k, l, kernel_row, kernel_col);
+ intermediate_sum += arr_mat_rm[index] * kernel_rm[kernel_index];
+ }
+ }
+ // Store the sum in the output matrix.
+ d_out_mat[index_to_row_major(i, j, row_eff, col_eff)] = intermediate_sum;
+ }
+ }
+}
+
+// __global__ void d_cuda_convolution(int* d_out_mat, int* arr_mat_rm, int* kernel_rm, int row_eff, int col_eff, int kernel_row, int kernel_col) {
+// // Calculate real row and column of input matrix.
+// int row = row_eff + kernel_row - 1;
+// int col = col_eff + kernel_col - 1;
+//
+// // Get i, and j from threadIdx
+// int tid = blockIdx.x * blockDim.x + threadIdx.x;
+// int i, j;
+// d_row_major_to_index(tid, row_eff, col_eff, &i, &j);
+//
+// // Calculate element in input matrix that is not on the boundary,
+// if (i < row_eff && j < col_eff) {
+// int intermediate_sum = 0;
+// for (int k = 0; k < kernel_row; k++) {
+// for (int l = 0; l < kernel_col; l++) {
+// int index = index_to_row_major(i + k, j + l, row, col);
+// int kernel_index = index_to_row_major(k, l, kernel_row, kernel_col);
+// intermediate_sum += arr_mat_rm[index] * kernel_rm[kernel_index];
+// }
+// }
+// d_out_mat[index_to_row_major(i, j, row_eff, col_eff)] = intermediate_sum;
+// }
+// }
+
// main() driver
@@ -249,10 +378,62 @@ int main() {
Matrix* arr_mat = (Matrix*)malloc(num_targets * sizeof(Matrix));
int arr_range[num_targets];
- // read each target matrix, compute their convolution matrices, and compute their data ranges
+ // Calculate variable for cuda computing.
+ int a = (target_row-kernel_row+1) * (target_col-kernel_col+1);
+ int b = 1024;
+ int block_size = a/b + (a % b != 0); // ceil(a/b)
+ int threads_per_block = 1024;
+ int row_eff = target_row - kernel_row + 1;
+ int col_eff = target_col - kernel_col + 1;
+
+ // Initialize host and device input and output matrixes.
+ int ** arr_mat_rm, **h_out_mat, ** d_out_mat, *kernel_rm;
+ // Allocate input matrix.
+ arr_mat_rm = (int**)malloc(sizeof(int*) * num_targets);
+ for (int i = 0; i < num_targets; i++) {
+ arr_mat_rm[i] = (int*)malloc(sizeof(int) * target_row * target_col);
+ }
+ // Allocate output matrix.
+ h_out_mat = (int**)malloc(sizeof(int*) * num_targets);
+ for (int i = 0; i < num_targets; i++) {
+ h_out_mat[i] = (int*)malloc(sizeof(int) * row_eff * col_eff);
+ }
+ // cudaMalloc((void**)&d_out_mat, sizeof(int*) * num_targets);
+ // for (int i = 0; i < num_targets; i++) {
+ // cudaMalloc(&h_out_mat[i], sizeof(int) * row_eff * col_eff);
+ // }
+ // cudaMemcpy(d_out_mat, h_out_mat, sizeof(int*) * num_targets, cudaMemcpyHostToDevice);
+
+ d_out_mat = (int**)malloc(sizeof(int*) * num_targets);
+ for (int i = 0; i < num_targets; i++) {
+ d_out_mat[i] = (int*)malloc(sizeof(int) * row_eff * col_eff);
+ }
+ kernel_rm = (int*)malloc(sizeof(int) * kernel_col * kernel_row);
+
+
+
+ // Store kernel in row major form.
+ kernel_rm = map_matrix(kernel.mat, kernel_row, kernel_col);
+
+ // read each target matrix, and get the row major matrix from.
for (int i = 0; i < num_targets; i++) {
arr_mat[i] = input_matrix(target_row, target_col);
- arr_mat[i] = convolution(&kernel, &arr_mat[i]);
+ arr_mat_rm[i] = map_matrix(arr_mat[i].mat, target_row, target_col);
+ // cuda_convolution<<<block_size, threads_per_block>>>(d_out_mat[i], arr_mat_rm[i], kernel_rm, target_row, target_col, kernel_row, kernel_col);
+ // cudaMemcpy(h_out_mat[i], d_out_mat[i], sizeof(int) * row_eff * col_eff, cudaMemcpyDeviceToHost);
+ cuda_convolution(d_out_mat[i], arr_mat_rm[i], kernel_rm, row_eff, col_eff, kernel_row, kernel_col);
+ arr_mat[i] = rm_to_matrix_object(d_out_mat[i], row_eff, col_eff);
+ }
+
+ // // Free cuda memory
+ // for (int i = 0; i < num_targets; i++) {
+ // cudaFree(h_out_mat[i]);
+ // }
+ // cudaFree(d_out_mat);
+
+ // For each target matrix, compute their convolution matrices, and compute their data ranges
+ for (int i = 0; i < num_targets; i++) {
+ // arr_mat[i] = convolution(&kernel, &arr_mat[i]);
arr_range[i] = get_matrix_datarange(&arr_mat[i]);
}