## 3 Our Proposed Approach
### GPU Kernel Optimization
However, CSE may degrade warp scheduling performance since CSE tends to use additional registers. We set four options according to the degree of CSE: low (L), medium-low (ML), medium-high (MH), and high (H). As the degree of CSE increases, more repeated expressions are removed by using more registers. Figure 6 shows a part of the source code of the L and H options. In the source code, the L option compute repeated expressions with fewer variable usage. However, the H option precomputes the repeated expressions just once, but it needs 15 more registers than the L option. Table 3 shows the total number of register variables substituting repeated expressions for each degree of CSE. Constant memory is a read-only memory that all threads can access. We can place the variable in the CUDA constant memory by declaring a variable with the keyword \"_constant_\". The constant memory variables originally reside in global memory, but they are cached automatically in a per-SM constant cache. Since they are not modified during kernel execution, the cache coherence issue does not occur. Furthermore, since the constant cache is optimized for broadcasting a value to massive threads, we can effectively reduce the register usage and memory bottleneck by the constant caching. Therefore, we set two options for constant caching: using constant caching (\"with CC\") and not using constant caching (\"without CC\"). When constant caching is used, we copy an array of 14 variables globally referenced by all threads to constant memory space, and the kernel uses the array directly. On the other hand, when constant caching is not used, the kernel accesses those variables in the global memory space. Figure 5: CSE implementation example. \begin{table}
\begin{tabular}{c c} \hline \hline
**CSE Option** & **The Number of Precomputed Subexpressions** \\ \hline L & 24 \\ ML & 28 \\ MH & 40 \\ H & 43 \\ \hline \hline \end{tabular}
\begin{tabular}{c c} \hline \hline
**Without CSE** & **With CSE** \\ \hline
1 & double AI(double x, double y) (...) \\
2 & double A2(double x, double y) (...) \\
3 & double A3(double x, double y) (...) \\
4 & 4 \\
5 & imt main() \\
6 & \{ \\
7 & // Calculate common expressions repeatedly \\
8 & double B1 = AI(x, y) + A2(x, y) + A3(x, y); \\
9 & double B2 = AI(x, y) + A2(x, y) + A3(x, y); \\
10 & // - \\
11 & 11 \\
12 & double B1 = tmpA1 + tmpA2 + tmpA3; \\
13 & double B2 = tmpA1 + tmpA2 + tmpA3; \\
14 & // - \\
15 & \} \\ \hline \hline \end{tabular}
\end{table} 
Table 3: CSE option description. The thread block size is an important factor for determining the number of active warps per SM because the CUDA restricts the maximum number of resident threads and thread blocks per SM. Since many active warps hide the latencies of stalled instructions, the thread block size can affect the performance of the kernel. Therefore, we set four options of the thread block size, i.e., 64, 128, 256, and 512, based on a rule of thumb that threads per block should be a multiple of warp size [35] and theoretical occupancy derived by the CUDA occupancy calculator [49]. All possible combinations from the options of the three optimization techniques are \(4\times 2\times 4=32\). We conducted experiments to evaluate the performance of our optimization approach. First, we compared 32 combinations using efficiency and occupancy to estimate the effects of options on the instruction efficiency and thread-level parallelism. Then, based on the computation time, we evaluated the significance of the instruction efficiency and thread-level parallelism on the parallel implementation of the Okada model and select the fastest one as the best CUDA optimization configuration.