diff --git a/makefile b/makefile
index 3d8f75e64f66c3cde25b7c3e6ad255ce577e51e8..6ac3bb473db80598c5e45ed57b7ecae23b089b56 100644
--- a/makefile
+++ b/makefile
@@ -1,29 +1,36 @@
-SOURCE_SERIAL=./src/serial.c
+SOURCE_SERIAL=./src/serial_real.c
 SOURCE_PARALLEL=./src/parallel.cu
 SOURCE_TEST=./src/test.c
+SOURCE_PARALLEL_COLLAB=./src/parallel_collab.cu
 
 EXEC_SERIAL=main-serial
 EXEC_PARALLEL=main-parallel
 EXEC_TEST=test
+EXEC_PARALLEL_COLLAB=parallel-collab
 
+COMPILER_PARALLEL_COLLAB=nvcc
 COMPILER_PARALLEL=nvcc
 COMPILER_SERIAL=gcc
 
 
 
 # Compile program.
+compile-parallel-collab:
+	${COMPILER_PARALLEL_COLLAB} -rdc=true -o ./bin/${EXEC_PARALLEL_COLLAB} ${SOURCE_PARALLEL_COLLAB}
 compile-serial:
 	${COMPILER_SERIAL} -o ./bin/${EXEC_SERIAL} ${SOURCE_SERIAL}
 compile-parallel:
 	${COMPILER_PARALLEL} -rdc=true -o ./bin/${EXEC_PARALLEL} ${SOURCE_PARALLEL} ./src/lib/matrix.cu ./src/lib/utils.cu ./src/lib/bitonic_sort.cu
 
 # Link program.
+link-parallel-collab : compile-parallel-collab
 link-serial: compile-serial
 link-parallel: compile-parallel
 
 # Compile and link.
 install-serial: compile-serial link-serial
 install-parallel: compile-parallel link-parallel
+install-parallel-collab: compile-parallel-collab link-parallel-collab
 
 # Clean bin folder.
 clean:
@@ -35,9 +42,14 @@ run-serial: install-serial
 run-serial-example: install-serial
 	./bin/${EXEC_SERIAL} < ./testcase/Example-TC
 run-parallel: install-parallel
-	./bin/${EXEC_PARALLEL
+	./bin/${EXEC_PARALLEL}
 run-parallel-example: install-parallel
 	./bin/${EXEC_PARALLEL} < ./testcase/Example-TC
+run-parallel-collab: install-parallel-collab
+	./bin/${EXEC_PARALLEL_COLLAB}
+run-parallel-collab-example: install-parallel-collab
+	./bin/${EXEC_PARALLEL_COLLAB} < ./testcase/Example-TC
+
 
 # Example testcase.
 TC-example: 
@@ -65,6 +77,17 @@ TC3-parallel: install-parallel
 TC4-parallel: install-parallel
 	(./bin/${EXEC_PARALLEL} < ./testcase/K01-03-TC4) > ./result/K01-03-TC4_parallel.txt
 
-TC-parallel: TC1-parallel TC2-parallel TC3-parallel TC4-parallel
+# Generic parallel testcase.
+TC1-parallel-collab: install-parallel-collab
+	(./bin/${EXEC_PARALLEL_COLLAB} < ./testcase/K01-03-TC1) > ./result/K01-03-TC1_parallel_collab.txt
+TC2-parallel-collab: install-parallel-collab
+	(./bin/${EXEC_PARALLEL_COLLAB} < ./testcase/K01-03-TC2) > ./result/K01-03-TC2_parallel_collab.txt
+TC3-parallel-collab: install-parallel-collab
+	(./bin/${EXEC_PARALLEL_COLLAB}  < ./testcase/K01-03-TC3) > ./result/K01-03-TC3_parallel_collab.txt
+TC4-parallel-collab: install-parallel-collab
+	(./bin/${EXEC_PARALLEL_COLLAB} < ./testcase/K01-03-TC4) > ./result/K01-03-TC4_parallel_collab.txt
 
+TC-parallel: TC1-parallel TC2-parallel TC3-parallel TC4-parallel
+TC-parallel-collab: TC1-parallel-collab TC2-parallel-collab TC3-parallel-collab TC4-parallel-collab
 TC-generic: TC-serial TC-parallel
+TC-generic-collab: TC-serial TC-parallel TC-parallel-collab
diff --git a/result/K01-03-TC1_parallel.txt b/result/K01-03-TC1_parallel.txt
index adfa022b810f199d3098dcaa2f4ec2e70267c0cf..3c886f1c9a184d1d7834ee879c979181d0acee66 100644
--- a/result/K01-03-TC1_parallel.txt
+++ b/result/K01-03-TC1_parallel.txt
@@ -3,4 +3,4 @@
 10114197
 10323010
 
-Runtime: 0.160294 s
+Runtime: 0.162468 s
diff --git a/result/K01-03-TC1_parallel_collab.txt b/result/K01-03-TC1_parallel_collab.txt
new file mode 100644
index 0000000000000000000000000000000000000000..2dca374dbec3557fff071cc00dcca8c03f35ae16
--- /dev/null
+++ b/result/K01-03-TC1_parallel_collab.txt
@@ -0,0 +1,6 @@
+8539213
+11916317
+10114197
+10323010
+
+Runtime: 0.163408 s
diff --git a/result/K01-03-TC1_serial.txt b/result/K01-03-TC1_serial.txt
index e828e833ca9c9337772da9a27546c56263a3d0fb..72c26b79ee8880472338f127d0e6385f0fe84481 100644
--- a/result/K01-03-TC1_serial.txt
+++ b/result/K01-03-TC1_serial.txt
@@ -3,4 +3,4 @@
 10114197
 10323010
 
-Runtime: 0.016049 s
+Runtime: 0.015577 s
diff --git a/result/K01-03-TC2_parallel.txt b/result/K01-03-TC2_parallel.txt
index 5bb2408e1415c74b18d449148f66d05b4ff4534c..af4e4b027a1b913436fe0cd10c5ba86b997bf818 100644
--- a/result/K01-03-TC2_parallel.txt
+++ b/result/K01-03-TC2_parallel.txt
@@ -3,4 +3,4 @@
 37739803
 38222937
 
-Runtime: 0.960191 s
+Runtime: 0.857004 s
diff --git a/result/K01-03-TC2_parallel_collab.txt b/result/K01-03-TC2_parallel_collab.txt
new file mode 100644
index 0000000000000000000000000000000000000000..760c0b3bc863cdd6edb1e569d3ccc1ff6e938105
--- /dev/null
+++ b/result/K01-03-TC2_parallel_collab.txt
@@ -0,0 +1,6 @@
+35064588
+46265294
+37739803
+38222937
+
+Runtime: 0.446493 s
diff --git a/result/K01-03-TC2_serial.txt b/result/K01-03-TC2_serial.txt
index c794bb80d78e0c28ec6fce6b2eb0d20767e84af8..b07814f47036817b4b1fb9bccd626985e0622703 100644
--- a/result/K01-03-TC2_serial.txt
+++ b/result/K01-03-TC2_serial.txt
@@ -3,4 +3,4 @@
 37739803
 38222937
 
-Runtime: 0.742690 s
+Runtime: 0.715588 s
diff --git a/result/K01-03-TC3_parallel.txt b/result/K01-03-TC3_parallel.txt
index c19d3246e0be37fe7e3a1d5112e3fd2be183d2ef..67a2084690716e97273ec09707988d14119e6d0d 100644
--- a/result/K01-03-TC3_parallel.txt
+++ b/result/K01-03-TC3_parallel.txt
@@ -3,4 +3,4 @@
 23198319
 23380111
 
-Runtime: 1.101716 s
+Runtime: 0.984349 s
diff --git a/result/K01-03-TC3_parallel_collab.txt b/result/K01-03-TC3_parallel_collab.txt
new file mode 100644
index 0000000000000000000000000000000000000000..ee9fbc695fb195d8d1de1e8e256a4eb5a2b56b5e
--- /dev/null
+++ b/result/K01-03-TC3_parallel_collab.txt
@@ -0,0 +1,6 @@
+18815130
+28707695
+23198319
+23380111
+
+Runtime: 0.722501 s
diff --git a/result/K01-03-TC3_serial.txt b/result/K01-03-TC3_serial.txt
index 8edba32113e491f17d87a3f1179d4cdb8be8df61..72325755f7032cdde2b4b971d2f8b566c2f9678c 100644
--- a/result/K01-03-TC3_serial.txt
+++ b/result/K01-03-TC3_serial.txt
@@ -3,4 +3,4 @@
 23198319
 23380111
 
-Runtime: 0.884338 s
+Runtime: 0.817598 s
diff --git a/result/K01-03-TC4_parallel.txt b/result/K01-03-TC4_parallel.txt
index 9ee11a5f49187c982513d57daff3b38720bf4ef1..3fa6f6a6ce39c1bfa26456f30b9e81c742467714 100644
--- a/result/K01-03-TC4_parallel.txt
+++ b/result/K01-03-TC4_parallel.txt
@@ -3,4 +3,4 @@
 51451884
 51774352
 
-Runtime: 10.189295 s
+Runtime: 9.561980 s
diff --git a/result/K01-03-TC4_parallel_collab.txt b/result/K01-03-TC4_parallel_collab.txt
new file mode 100644
index 0000000000000000000000000000000000000000..b6c4758c2d5dc6d887a68228628bb237ae3ac078
--- /dev/null
+++ b/result/K01-03-TC4_parallel_collab.txt
@@ -0,0 +1,6 @@
+41250811
+71841136
+51451884
+51774352
+
+Runtime: 4.290069 s
diff --git a/result/K01-03-TC4_serial.txt b/result/K01-03-TC4_serial.txt
index d947b0b2f67a0d60574286b559ffb2c4ad8787a4..5740df2da2fa1e21ace34b7f63d2a29a9dddde29 100644
--- a/result/K01-03-TC4_serial.txt
+++ b/result/K01-03-TC4_serial.txt
@@ -3,4 +3,4 @@
 51451884
 51774352
 
-Runtime: 9.711713 s
+Runtime: 8.956075 s
diff --git a/src/lib/bitonic_sort.cu b/src/lib/bitonic_sort.cu
index 7b8c1816c29580c4435b1a3bfccf7e373dd731ad..b6d6d82c4842a0ab81ac13e0d3a16c316d77f4a3 100644
--- a/src/lib/bitonic_sort.cu
+++ b/src/lib/bitonic_sort.cu
@@ -126,6 +126,7 @@ __global__ void bitonic_sort_step(int *d_arr, int i, int j)
 
   // The thread index.
   idx = threadIdx.x + blockDim.x * blockIdx.x;
+
   // The thread index of the patner.
   patner = idx ^ j;
 
diff --git a/src/parallel.cu b/src/parallel.cu
index f43975e01362a5e1c4c54191dcfaa013aa432c1c..d6d0bb217003fe6d0798645722ca760a2e61db1c 100644
--- a/src/parallel.cu
+++ b/src/parallel.cu
@@ -9,7 +9,6 @@
 #include "./lib/utils.cuh"
 
 
-// main() driver
 int main() {
 	// Time.
 	clock_t t;
diff --git a/src/parallel_collab.cu b/src/parallel_collab.cu
index f214d43ebbe9e28218fc3992e82163c6b9b4e981..5d6470638f34414498fc42190687f246aa554870 100644
--- a/src/parallel_collab.cu
+++ b/src/parallel_collab.cu
@@ -1,4 +1,4 @@
-// parallel.cu
+// parallel_collab.cu
 
 #include <stdio.h>
 #include <stdlib.h>
@@ -137,58 +137,145 @@ Matrix convolution(Matrix *kernel, Matrix *target) {
 }
 
 
-/*
- * 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];
-	}
+/**
+ * Swap the values of two elements in an array.
+ * 
+ * @param d_arr - the array.
+ * @param i - the index of the first element.
+ * @param j - the index of the second element.
+ */
+__device__ void swap(int *d_arr, int i, int j)
+{
+  int temp = d_arr[i];
+  d_arr[i] = d_arr[j];
+  d_arr[j] = temp;
+}
 
-	for (int i = 0; i < n_right; i++) {
-		arr_right[i] = n[i + mid + 1];
-	}
+/**
+ * Pad the remaining empty value in the array with the maximum value.
+ * 
+ * @param arr - the array.
+ * @param length - the number of elements inside the array.
+ * @param length_buffer - the length of the buffer.
+ */
+void padding_array(int *arr, int length, int length_buffer)
+{
+  for (int i = length; i < length_buffer; i++) {
+    arr[i] = INT_MAX;
+  }
+}
 
-	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++;
-	}
+/**
+ * Copy the array from src to dest and pad the remaining empty value with the maximum value.
+ * 
+ * @param dest - the destination array.
+ * @param src - the source array.
+ * @param length - the number of elements inside the src array.
+ * @param length_buffer - the length of the buffer.
+ */
+void copy_padding(int *dest, int *src, int length, int length_buffer)
+{
+  for (int i = 0; i < length; i++) {
+    dest[i] = src[i];
+  }
+  padding_array(dest, length, length_buffer);
+}
 
-	while (iter_left < n_left)  {
-		n[iter_merged++] = arr_left[iter_left++];
-	}
-	while (iter_right < n_right) {
-		n[iter_merged++] = arr_right[iter_right++];
-	} 
+/**
+ * Get the number of minimum blocks as a power of 2.
+ * 
+ * @param length - the number of elements inside the array.
+ * @param num_threads - the number of threads.
+ * @return int - the number of minimum blocks needed.
+ */
+int minimum_blocks(int length, int num_threads) {
+  int num_blocks = 1;
+
+  // Increase the number of blocks 2 times, so it will be a power of 2.
+  while (num_blocks * num_threads < length) {
+    num_blocks *= 2;
+  }
+
+  return num_blocks;
 }
 
+/**
+ * Do the bitnoic sort step by step.
+ * 
+ * @param d_values - array in the device to be sorted.
+ * @param i - major step index.
+ * @param j - minor step index.
+ */
+__global__ void bitonic_sort_step(int *d_arr, int i, int j)
+{
+  // The array index and its patner.
+  int idx, patner;
+
+  // The thread index.
+  idx = threadIdx.x + blockDim.x * blockIdx.x;
+
+  // The thread index of the patner.
+  patner = idx ^ j;
+
+  // Sort the array by threads with the lowest idx.
+  if (idx < patner) {
+    if ((idx & i) == 0) {
+      // Sort ascending.
+      if (d_arr[idx] > d_arr[patner]) {
+        swap(d_arr, idx, patner);
+      }
+    }
+    if ((idx & i) != 0) {
+      // Sort descending.
+      if (d_arr[idx] < d_arr[patner]) {
+        swap(d_arr, idx, patner);
+      }
+    }
+  }
+}
 
-/* 
- * Procedure merge_sort
+/**
+ * Perform a bitonic sort on the array.
  *
- * 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;
+ * @param h_arr The host array to sort.
+ * @param length The length of the array.
+ */
+void bitonic_sort(int *h_arr, int length)
+{
+  // Initialize the constants variable.
+  const int threads = 1024;
+  const int blocks = minimum_blocks(length, threads);
+  const int buffer_length = threads * blocks;
+
+  // Initialize the memory size of the array.
+  size_t size = length * sizeof(int);
+  size_t buffer_size = buffer_length * sizeof(int);
+
+  // Create the buffer array and pad with maximum value of Int.
+  int *h_buffer = (int *)malloc(buffer_size);
+  copy_padding(h_buffer, h_arr, length, buffer_length);
+
+  // Allocate and copy array into device memory.
+  int *d_arr;
+  cudaMalloc((void**) &d_arr, buffer_size);
+  cudaMemcpy(d_arr, h_buffer, buffer_size, cudaMemcpyHostToDevice);
+
+  // Sort the array using bitonic_sort_step.
+  int i, j;
+  // The major step.
+  for (i = 2; i <= buffer_length; i *= 2) {
+    // The minor step.
+    for (j = i / 2; j > 0; j = j / 2) {
+      bitonic_sort_step<<<blocks, threads>>>(d_arr, i, j);
+    }
+  }
 
-		merge_sort(n, left, mid);
-		merge_sort(n, mid + 1, right);
+  // Copy the values back to the host.
+  cudaMemcpy(h_arr, d_arr, size, cudaMemcpyDeviceToHost);
 
-		merge_array(n, left, mid, right);
-	}	
+  // Free device memory.
+  cudaFree(d_arr);
+  free(h_buffer);
 }
  
 
@@ -240,7 +327,7 @@ int index_to_row_major(int row, int col, int row_eff, int col_eff) {
 }
 
 __device__ int d_index_to_row_major(int row, int col, int row_eff, int col_eff) {
-		return row * col_eff + col;
+ 		return row * col_eff + col;
 }
 
 /**
@@ -273,6 +360,22 @@ int* map_matrix(int mat[][100], int row, int col) {
 	return map;
 }
 
+/*
+ * Function map_matrix_extended
+ * 
+ * Returns a row major matrix of the input matrix.
+ **/
+int* map_matrix_extended(int** mat, 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
  * 
@@ -305,68 +408,117 @@ Matrix rm_to_matrix_object(int* map, int row, int col) {
 	return mat;
 }
 
+/**
+ * Function rm_to_list_matrix_object
+ * 
+ * Return List of Matrix Struct of row major matrix
+ */
+Matrix* rm_to_list_matrix_object(int* map, int row, int col, int row_inner, int col_inner) {
+	Matrix* mat = (Matrix*) malloc(sizeof(Matrix) * row);
+	for (int i = 0; i < row; i++) {
+		init_matrix(&mat[i], row_inner, col_inner);
+		int pad = i * col;
+		for (int j = 0; j < row_inner; j++) {
+			for (int k = 0; k < col_inner; k++) {
+				int index = index_to_row_major(j, k, row_inner, col_inner) + pad;
+				mat[i].mat[j][k] = map[index];
+			}
+		}
+	}
+	return mat;
+}
+/**
+ * Function list_matrix_object_to_rm
+ * 
+ * Return row major matrix of list of Matrix struct
+ */
+int* list_matrix_object_to_rm(Matrix* mat, int num_matrix,  int row_inner, int col_inner) {
+	int* map = (int*) malloc(sizeof(int) * num_matrix * row_inner * col_inner);
+	for (int i = 0; i < num_matrix; i++) {
+		int pad = i * row_inner * col_inner;
+		for (int j = 0; j < row_inner; j++) {
+			for (int k = 0; k < col_inner; k++) {
+				int index = index_to_row_major(j, k, row_inner, col_inner) + pad;
+				map[index] = mat[i].mat[j][k];
+			}
+		}
+	}
+	
+	return map;
+}
+
 /**
  * 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) {
+void cuda_convolution(int* out_mat_rm, int* arr_mat_rm, int* kernel_rm, int row_eff, int col_eff, int kernel_row, int kernel_col, int curr_mat) {
 	// Calculate real row and column of input matrix.
 	int row = row_eff + kernel_row - 1;
 	int col = col_eff + kernel_col - 1;
 	
+	// Calculate padding target and output matrix.
+	int pad = curr_mat * row * col;
+	int pad_out = curr_mat * row_eff * col_eff;
+
 	// For each element in input matrix that is not on the boundary,
-	for (int i = 0; i < row_eff; i++) {
+	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 index = index_to_row_major(i + k, j + l, row, col) + pad;
 					int kernel_index = index_to_row_major(k, l, kernel_row, kernel_col);
 					intermediate_sum += arr_mat_rm[index] * kernel_rm[kernel_index];
+					// Print all i,j,k,l
+					// printf("i:%d, j:%d, k:%d, l:%d\n", i, j, k, l);
+
 				}
 			}
 			// Store the sum in the output matrix.
-			d_out_mat[index_to_row_major(i, j, row_eff, col_eff)] = intermediate_sum;
+			out_mat_rm[index_to_row_major(i, j, row_eff, col_eff) + pad_out ] = 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) {
-  printf("aaaa\n");
-	// Calculate real row and column of input matrix.
- 	int row = row_eff + kernel_row - 1;
- 	int col = col_eff + kernel_col - 1;
+  // Calculate real row and column of input matrix.
+	int row = row_eff + kernel_row - 1;
+	int col = col_eff + kernel_col - 1;
+
+  // Determine current matrix from block;
+	int curr_mat = blockIdx.y;
+  // Calculate padding target and output matrix.
+ 	int pad = curr_mat * row * col;
+ 	int pad_out = curr_mat * row_eff * col_eff;
  
   // 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);
-	printf("bbbb\n");
-
-  // Calculate element in input matrix that is not on the boundary
+ 
+  // 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 = d_index_to_row_major(i + k, j + l, row, col);
+ 				int index = d_index_to_row_major(i + k, j + l, row, col) + pad;
  				int kernel_index = d_index_to_row_major(k, l, kernel_row, kernel_col);
  				intermediate_sum += arr_mat_rm[index] * kernel_rm[kernel_index];
  			}
  		}
- 		d_out_mat[d_index_to_row_major(i, j, row_eff, col_eff)] = intermediate_sum;
+ 		d_out_mat[d_index_to_row_major(i, j, row_eff, col_eff) + pad_out] = intermediate_sum;
  	}
 }
 
 
-
 // main() driver
 int main() {
 	// Time.
 	clock_t t;
-  t = clock();
+    t = clock();
 
 	int kernel_row, kernel_col, target_row, target_col, num_targets;
 	
@@ -388,68 +540,85 @@ int main() {
 	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);
+
+	// Allocate variable. 
+	// rm means row-major. It's indicate matrix are in row-major order.
+	// Variable declaration.
+	int * arr_mat_rm=0, * d_arr_mat_rm=0, *out_mat_rm=0, *d_out_mat_rm=0, *kernel_rm=0, *d_kernel_rm=0;
+	int size_arr_mat, size_out_mat, size_kernel;
+	cudaError err;
+
+	// Allocate input matrix in device and host.
+	size_arr_mat = num_targets * target_row * target_col;
+	arr_mat_rm = (int*)malloc(sizeof(int) * size_arr_mat);
+	cudaMalloc((void **)&d_arr_mat_rm, sizeof(int) * size_arr_mat);
+	if (arr_mat_rm == 0 | d_arr_mat_rm == 0) {
+		printf("Error: Memory allocation failed for arr_mat.\n");
+	 	return 1;
 	}
-	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);
-		}
-	cudaError err = cudaMemcpy(d_out_mat, h_out_mat, sizeof(int*) * num_targets, cudaMemcpyHostToDevice);
-	if(err!=cudaSuccess) {
-    printf("CUDA error copying to Host: %s\n", cudaGetErrorString(err));
-  }
 
-	// 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);
-	
+	// Allocate output matrix in device and host.
+	size_out_mat = num_targets * row_eff * col_eff;
+	out_mat_rm = (int*)malloc(sizeof(int) * size_out_mat);
+	cudaMalloc((void **)&d_out_mat_rm, sizeof(int) * size_out_mat);
+	if (out_mat_rm == 0 | d_out_mat_rm == 0) {
+	 	printf("Error: Memory allocation failed for out_mat.\n");
+	 	return 1;
+	}
+	cudaMemset(d_out_mat_rm, 0, sizeof(int) * size_out_mat);
 
 
-	// Store kernel in row major form.
+	// Allocate kernel matrix in host.
+	size_kernel = kernel_row * kernel_col;
+	kernel_rm = (int*)malloc(sizeof(int) * size_kernel);
+	// Store kernel in row major form and allocate kernel for device.
 	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++) {			
-		printf("a\n");
+	cudaMalloc((void **)&d_kernel_rm, sizeof(int) * size_kernel);
+	if (kernel_rm == 0 | d_kernel_rm == 0) {
+	 	printf("Error: Memory allocation failed for kernel.\n");
+	 	return 1;
+	}
+	err = cudaMemcpy(d_kernel_rm, kernel_rm, sizeof(int) * size_kernel, cudaMemcpyHostToDevice);
+	 if (err != cudaSuccess) {
+		 	printf("Error copy host to device: %s\n", cudaGetErrorString(err));
+			return 1;
+	}
+	
+	
+	// Read each target matrix.
+	for (int i = 0; i < num_targets; i++) {
 		arr_mat[i] = input_matrix(target_row, target_col);
-		printf("b\n");
-		arr_mat_rm[i] = map_matrix(arr_mat[i].mat, target_row, target_col);
-		printf("c\n");
-		d_cuda_convolution<<<block_size, threads_per_block>>>(d_out_mat[i], arr_mat_rm[i], kernel_rm, row_eff, col_eff, kernel_row, kernel_col);
-		printf("d\n");
-		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);
-		printf("e\n");
-		arr_mat[i] = rm_to_matrix_object(d_out_mat[i], row_eff, col_eff);
-	  printf("f\n");
+	}
+	// Store each target matrix in row major form and allocate target matrix for device.
+	arr_mat_rm = list_matrix_object_to_rm(arr_mat, num_targets, target_row, target_col);
+	err = cudaMemcpy(d_arr_mat_rm, arr_mat_rm, sizeof(int) * size_arr_mat, cudaMemcpyHostToDevice);
+	if (err != cudaSuccess) {
+		printf("Error copy host to device: %s\n", cudaGetErrorString(err));
+		return 1;
 	}
 
-	// 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.
+	dim3 grid, block;
+	grid.x = block_size;
+	grid.y = num_targets;
+	block.x = threads_per_block;
+	d_cuda_convolution<<<grid, block>>>(d_out_mat_rm, d_arr_mat_rm, d_kernel_rm, row_eff, col_eff, kernel_row, kernel_col);
+	err = cudaMemcpy(out_mat_rm, d_out_mat_rm, sizeof(int) * size_out_mat, cudaMemcpyDeviceToHost);
+	if (err != cudaSuccess) {
+	 	printf("Error copy device to host: %s\n", cudaGetErrorString(err));
+	}
+	// for (int i = 0; i < num_targets; i++){
+		// cuda_convolution(out_mat_rm, arr_mat_rm, kernel_rm, row_eff, col_eff, kernel_row, kernel_col, i);
+	// }
+	arr_mat = rm_to_list_matrix_object(out_mat_rm, num_targets, row_eff*col_eff, row_eff, col_eff);
 
-	// For each target matrix, compute their convolution matrices, and compute their data ranges
+	// For each target matrix, 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]); 
 	}
 
 	// sort the data range array
-	merge_sort(arr_range, 0, num_targets - 1);
+	bitonic_sort(arr_range, num_targets);
 	
 	int median = get_median(arr_range, num_targets);	
 	int floored_mean = get_floored_mean(arr_range, num_targets); 
@@ -463,7 +632,11 @@ int main() {
 
 	// Print execution time in seconds.
 	t = clock() - t;
-	printf("\nRuntime: %f s\n", ((float)t) / CLOCKS_PER_SEC);
+	printf("\nRuntime: %f s\n", ((float)t) / CLOCKS_PER_SEC);	
 	
+	// Free cuda Memory.
+	cudaFree(d_arr_mat_rm);
+	cudaFree(d_out_mat_rm);
+	cudaFree(d_kernel_rm);
 	return 0;
-}
\ No newline at end of file
+}
diff --git a/src/serial _real.c b/src/serial_real.c
similarity index 95%
rename from src/serial _real.c
rename to src/serial_real.c
index ab02e0f5a4eb458acf82d82df7410930d60ba64f..7a641e5c891cbd5550fef486d2e0e47d9f13d45f 100644
--- a/src/serial _real.c	
+++ b/src/serial_real.c
@@ -1,277 +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;
-}
+// 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;
+}