# Changes

## Group 6

, 22:00, 7 April 2019
m
Assignment 3 - Optimize
# [mailto:xhuang110@myseneca.ca?subject=gpu610 Xiaowei Huang]
# [mailto:yyuan34@myseneca.ca?subject=gpu610 Yihang Yuan]
# [mailto:zzhou33@myseneca.ca?subject=dps915 Zhijian Zhou]
[mailto:xhuang110@myseneca.ca,yyuan34@myseneca.ca,zzhou33@myseneca.ca?subject=dps901-gpu610 Email All]
== Progress ==
==== The Monte Carlo Simulation (PI Calculation) ====
Subject: The Monte Carlo Simulation (PI Calculation)
Got the code from here:
https://rosettacode.org/wiki/Monte_Carlo_methods#C.2B.2B
It uses random sampling to define constraints on the value and then makes a sort of "best guess."
{| class="wikitable mw-collapsible mw-collapsed"! Source Code|-|<pre>#include<iostream>#include<fstream>#include<math.h>#include<stdlib.h>#include<time.h> using namespace std; void calculatePI(int n, float* h_a) { float x, y; int hit; srand(time(NULL)); for (int j = 0; j < n; j++) { hit = 0; x = 0; y = 0; for (int i = 0; i < n; i++) { x = float(rand()) / float(RAND_MAX); y = float(rand()) / float(RAND_MAX); if (y <= sqrt(1 -(x * x))) { hit += 1; } } h_a[j] = 4 * float(hit) / float(n); }} int main(int argc, char* argv[]) {  if (argc != 2) { std::cerr << argv[0] << ": invalid number of arguments\n"; std::cerr << "Usage: " << argv[0] << " size_of_matrices\n"; return 1; } int n = std::atoi(argv[1]); // scale float* cpu_a; cpu_a = new float[n];   calculatePI(n, cpu_a);  ofstream h_file; h_file.open("h_result.txt"); float cpuSum = 0.0f; for (int i = 0; i < n; i++) { cpuSum += cpu_a[i]; h_file << "Host: " << cpu_a[i] << endl; } cpuSum = cpuSum / (float)n; cout << "CPU Result: " << cpuSum << endl; h_file.close();}  </pre>|} As this algorithm is based on random sampling, so there is only one function that does all the work.Flat profile:<pre>Each sample counts as 0.01 seconds. % cumulative self self total time seconds seconds calls Ts/call Ts/call name 101.05 0.02 0.02 calculatePI(int, float*) 0.00 0.02 0.00 1 0.00 0.00 _GLOBAL__sub_I__Z11calculatePIiPf</pre> Call graph<pre>granularity: each sample hit covers 2 byte(s) for 0.47% of 2.11 seconds index % time self children called name <spontaneous>[1] 100.0 2.11 0.00 calculatePI(int, float*) [1]----------------------------------------------- 0.00 0.00 1/1 __libc_csu_init [17][9] 0.0 0.00 0.00 1 _GLOBAL__sub_I__Z11calculatePIiPf [9]-----------------------------------------------Index by function name  [9] _GLOBAL__sub_I__Z11calculatePIiPf (a1.cpp) [1] calculatePI(int, float*)</pre> '''Results for different scale of calculation'''
[[File:Yihang.JPG]]
[[File:yihang_p2_profile.png]]
----
=== Assignment 2 - Parallelize ===
In conclusion, by parallelized the serial version of the algorithm, we see an immediate improvement of performance.
----
=== Assignment 3 - Optimize ===
The main kernel optimization did on 2 parts

The new '''third kernel''' is to sum up the PI results which generated in each block, by applying serial reduction algorithm which could reduce its Big-O classification is O(log n). Beside, it also sets the limitation to check whether the calculation is out of the range. As it can stop the sum-up to those useless thread value. And the total value of that block is passed to the first element of that block.
<pre>
// kernel 3
// the third kernel sum the result in each block
__global__ void sumPi(float* d_a, float*d_b, const int n) {
int i = blockIdx.x * blockDim.x + threadIdx.x;
__shared__ float s[ntpb];
s[t] = d_a[i];

// sum the data in shared memory
for (int stride = 1; stride < blockDim.x; stride <<= 1) {
if ((t % (2 * stride) == 0) && (i + stride < n)) {
s[t] += s[t + stride];
}
}

// store the sum in d_b;
if (t == 0) {
d_b[blockIdx.x] = s[0];
}
}
</pre>

To reduce the data transfer from device to host, the new '''forth kernel''' is created. It uses to sum up all blocks results from passing back to the host, which also applies the reduction algorithm.
<pre>
// kernel 4
// the forth kernel sum the result of all blocks
__global__ void accumulate(float* c, const int nblocks) {
// store the elements of c[] in shared memory
int i = blockIdx.x * blockDim.x + threadIdx.x;
__shared__ float s[ntpb];
s[t] = c[i];

// sum the data in shared memory
for (int stride = 1; stride < blockDim.x; stride <<= 1) {
if ((t % (2 * stride) == 0) && (i + stride < nblocks)) {
s[t] += s[t + stride];
}
}

// store the sum in c[0]
if (t == 0) {
c[blockIdx.x] = s[0];
}
}
</pre>

Even though the runtime is closed after implementing the optimization, the kernels has a heavy load than before. It also finishes the accumulation of all randomly generated PI. I consider this optimization performs well.

[[File:Optimized PI calculation and accumulation.jpg]]
57
edits