14 This paper presents a Numba-based solver for the 1D Heat Equation, seamlessly 15 blending Python’s readability with Numba’s dynamic Just-In-Time (JIT) compi- 16 lation. The explicit method exhibits a notable runtime reduction from 8.324 s to 17 4.035 s, while the implicit method sees a more pronounced improvement, decreas- 18 ing from 9.970 s to 1.195 s. Statistical tests confirm the statistical significance of 19 these efficiency gains. 20 Future research directions include extending the solver to multidimensional 21 heat equations, exploring advanced parallelization techniques, and implementing 22 dynamic parameter optimization strategies. Collaboration with domain experts 23 for real-world applications is also envisioned to validate the solver’s performance 24 and impact. 25 In summary, the symbiosis of Python and Numba in crafting an optimized 1D 26 Heat Equation solver marks a pivotal advancement in efficient numerical solu- 27 tions. This research holds promise for diverse scientific applications, ushering in 28 a new era of computational efficiency. 29 Keywords: Computational Efficiency, Finite Difference Methods, Numerical PDE 30 Solver, Numba Optimization 31 1 1 Introduction 32 The relentless pursuit of efficient numerical solutions to partial differential equations 33 (PDEs) has been a catalyzing force propelling significant strides in the field of scien- 34 tific computing [1–3]. Amidst the pantheon of these equations, the 1D Heat Equation 35 emerges as a foundational model, wielding far-reaching applications in diverse scien- 36 tific realms, including physics and engineering. The precise simulation of temporal 37 heat evolution within a one-dimensional space is of paramount importance, providing 38 profound insights into an array of phenomena, ranging from the conductive intricacies 39 of materials to the nuanced thermal processes transpiring within electronic devices [e. 40 g., 4–8]. 41 In response to this imperative, we navigate the intricate landscape of numerical 42 computation, spotlighting the critical need for a bespoke solver tailored specifically 43 for the 1D Heat Equation. Our focus converges on harnessing the transformative 44 capabilities of Numba, a dynamic Just-In-Time (JIT) compiler meticulously designed 45 for the Python programming language [9]. Python, with its unparalleled readability, 46 versatility, and expansive library ecosystem, has rightfully claimed its dominance in 47 scientific computing [10, 11]. However, the inherent interpretive nature of Python poses 48 challenges to computational speed. Enter Numba—an ingenious JIT compiler that 49 dynamically translates Python functions into machine code during runtime, thus con- 50 ferring a substantial performance boost. It is the symbiosis between Python’s inherent 51 elegance and Numba’s computational alacrity that forms the nucleus of our mission: to 52 forge a solver that not only encapsulates the aesthetic appeal of Python programming 53 but also adeptly exploits the efficiency afforded by Numba’s JIT compilation. 54 The numerical solution of PDEs is an intricate dance with algorithms, where 55 optimization is the linchpin, striking a harmonious balance between accuracy and com- 56 putational efficiency. Crafting an optimized solver for the 1D Heat Equation mandates 57 a thoughtful consideration of numerical stability, precision, and algorithmic efficiency. 58 This paper embarks on a comprehensive expedition into the labyrinth of challenges 59 encountered during the optimization process, meticulously unfolding the narrative of 60 how Numba’s capabilities are judiciously harnessed to surmount these challenges and 61 elevate the solver’s overall performance. 62 A pivotal facet of our contribution lies in the granular presentation—a step-by- 63 step expos´ e of our Numba-based solver for the 1D Heat Equation. Here, we unveil the 64 intricacies of implementation, expound upon optimization strategies, and subject the 65 solver to rigorous performance evaluations. Beyond its practical utility as a robust 66 tool for researchers and practitioners immersed in heat transfer simulations, our work 67 extends its tendrils into the broader discourse on optimizing numerical methods in 68 Python. By illuminating the subtleties of our approach, we aspire to contribute to 69 the evolving understanding of how Numba can be judiciously employed to expedite 70 scientific computations, thereby forging a path towards more efficient and scalable 71 solutions to PDEs within the Python programming paradigm. 72 2 2 Methods 73 2.1 Numerical Experiments 74 The general form of the 3D heat equation was considered: 75 ∂T ∂t=α∇2T (1) , where Trepresents the temperature distribution in a three-dimensional space as a 76 function of both spatial coordinates ( x, y, z ) and time ( t),tdenotes time, αrepresents 77 the thermal diffusivity of the material, ∇2is the Laplacian operator, defined as ∇2T≡ 78 ∂2T ∂x2+∂2T ∂y2+∂2T ∂z2. 79 The case was considered where the temperature distribution Twas assumed to be 80 solely a function of one spatial coordinate, say x. Mathematically, T≡T(x, t). This 81 assumption simplified the Laplacian operator to only the second spatial derivative: 82 ∇2T≡∂2T ∂x2(2) The simplified Laplacian for 1D was substituted into the 3D heat equation: 83 ∂T ∂t=α∂2T ∂x2(3) The derived equation represents the 1D Heat equation. It describes the evolution of 84 the temperature distribution in a material along a single spatial dimension ( x) as a 85 function of time ( t) due to heat conduction. This equation is classified as a parabolic 86 PDE. 87 In this study, equation (3) was numerically solved using both explicit and implicit 88 methods of the finite difference approach, a fundamental technique widely employed 89 for approximating numerical solutions of PDEs [12]. Additionally, the results obtained 90 through these methods were analyzed to gain insights into the behavior of the system. 91 In this study, we employed a discretization approach to model the behavior of 92 temperature within a spatial domain. The spatial discretization involves partitioning 93 the domain into grid points, each separated by a uniform grid spacing ∆ x. Specifically, 94 we denote the temperature at each spatial point xiand time tnasTn i. 95 For the temporal discretization using the explicit scheme, we utilized the forward 96 difference for time. The update formula is expressed as: 97 Tn+1 i=Tn i+C(Tn i+1−2Tn i+Tn i−1) (4) , where C=α·∆t ∆x2, representing the Courant–Friedrichs–Lewy (CFL) number [13, 14]. 98 The stability of the explicit scheme is governed by the CFL condition, ensuring that 99 C≤1 2for numerical stability. 100 On the other hand, the implicit scheme employs the backward difference for time, 101 leading to the implicit update formula: 102 Tn+1 i=Tn i+C(Tn+1 i+1−2Tn+1 i+Tn+1 i−1) (5) 3 The implicit scheme, while unconditionally stable, involves solving a system of 103 equations at each time step, introducing additional computational cost. To express 104 the implicit scheme in a matrix form, we defined a tridiagonal matrix Arelated to the 105 discretization of the second spatial derivative. The matrix equation becomes: 106 ATn+1=Tn(6) , where Tn+1is the vector of temperatures at the next time step, Tnis the vector at 107 the current time step, and Ais the tridiagonal matrix. The tridiagonal matrix Awas 108 constructed as follows: 109 A= 1 + 2 C−C 0··· 0 −C1 + 2 C−C··· 0 0 −C1 + 2 C...... ............−C 0 0 ··· − C1 + 2 C (7) Each element Aijin the matrix represents the coefficients associated with the 110 temperatures Tn+1 iin the implicit update formula. This comprehensive methodology 111 enables the numerical simulation of temperature evolution within the defined spatial 112 domain. 113 In this study, the automated matrix calculation process utilized the NumPy library 114 [15]. Dirichlet boundary conditions were applied to both sides, with Tleftset to 20◦C 115 andTrightto 100◦C for demonstration purposes (Figure 1). The material chosen was 116 Plexiglass, characterized by a thermal diffusivity of 0.2 m2/s [16–18]. Users retain the 117 flexibility to adjust these values as needed. 118 Fig. 1 : Schematic depicting heat transfer in a 1D cross-section, featuring Dirichlet boundary conditions at both ends. In this demonstration, the spatial range extends from 0 to 1 m, with a spatial step 119 size of 0.01 m. Employing the explicit method, a time step size of 0.00025 s was used 120 within the temporal range from 0 seconds to 1 s, adhering to the upper limit of CFL 121 4 stability criterion. Users are free to select their preferred time step size for the implicit 122 method. 123 In the Numba-based simulation, we applied the 124 @jit(nopython=True, parallel=True) decorator to optimize the numerical calcu- 125 lation scheme function. This decorator enabled the translation of the code directly 126 into machine language via LLVM [19], facilitating parallelization of the calculation 127 process. However, due to the current limitations of Numba, we chose not to parallelize 128 the image plotting process, particularly with dynamic functions from the Matplotlib 129 library [20]. 130 2.2 Statistical Analysis of the Total Runtime 131 In our experimental design, encompassing both control and test experiments, the for- 132 mer executed the standard numerical algorithm without Numba optimizations, while 133 the latter implemented the Numba-optimized version. Each set of experiments under- 134 went 100 iterations to bolster statistical robustness. To quantify the impact of Numba 135 optimization on computational efficiency, we employed the hyperfine command-line 136 benchmarking tool [21]. This tool facilitated precise measurement of execution times 137 for each iteration, offering a comprehensive understanding of Numba’s discernible 138 influence on overall computational efficiency. 139 For our statistical analysis, we chose nonparametric tests to accommodate the 140 nature of the data and to ensure robustness against potential deviations from normality 141 assumptions. 142 To compare execution times between control and test experiments, we utilized 143 the Mann-Whitney U test [22], a nonparametric alternative suitable for scenarios 144 where normality assumptions may not be satisfied. The data from both samples were 145 combined and ranked from smallest to largest. Let U1be the sum of ranks for the 146 total runtime for control group, and U2be the sum of ranks for the total runtime of 147 the test group. The U statistic, given by U= min( U1, U2), provided a measure of the 148 difference between the groups. The expected value of U ( E(U)) was calculated as: 149 E(U) =n1n2 2(8) This assessed the central tendency under the null hypothesis. The standard deviation 150 of U ( SD(U)), given by: 151 SD(U) =r n1n2(n1+n2+ 1) 12(9) , provided a measure of dispersion. The z-score, calculated as: 152 z=U−E(U) SD(U)(10) , was used to derive the p-value from a standard normal distribution. A significance 153 level of 0.01 was chosen for the test, and if p <0.01, the null hypothesis was rejected. 154 For the paired nature of our experiments, we employed the Wilcoxon signed-rank 155 test [23], a nonparametric test suitable for comparing paired samples when normality 156 assumptions may not hold. The absolute differences between paired observations were 157 5 ranked, and positive and negative ranks were summed separately. Let W+be the sum 158 of ranks for positive differences, and W−be the sum of ranks for negative differences. 159 The test statistic W, calculated as: 160 W= min( W+, W−) (11) , represented the directionality of differences. The expected value of W ( E(W)), given 161 by: 162 E(W) =n(n+ 1) 4(12) , where nis the number of pairs, assessed the central tendency under the null 163 hypothesis. The standard deviation of W ( SD(W)), calculated as: 164 SD(W) =r n(n+ 1)(2 n+ 1) 24(13) , provided a measure of dispersion. The z-score, given by: 165 z=W−E(W) SD(W)(14) , was used to determine the p-value from a standard normal distribution. Adopting a 166 significance level of 0.01, the null hypothesis was rejected if p <0.01. 167 The utilization of these nonparametric tests and their derivations played a pivotal 168 role in establishing dependable statistical inferences during our evaluation of Numba 169 optimization’s effects. This ensured the resilience and credibility of our findings within 170 practical, real-world testing scenarios. Widely employed in the analysis of numerical 171 model outputs across diverse scientific domains [e. g., 24–28], both of these statistical 172 test methods were automatically implemented through the SciPy library [29]. 173 3 Results and Discussion 174 In Figure 2, the temporal evolution of temperature distribution is depicted utilizing 175 an explicit scheme with a fixed time step of ∆ t= 0.00025 s. This graph illustrates 176 the intricate changes in temperature over time, providing a detailed insight into the 177 dynamic behavior of the computational model within the specified time range from 178 0 to 10. The explicit scheme showcases a high level of accuracy and responsiveness, 179 capturing fine-grained temporal dynamics inherent in the 1D Heat Equation. The 180 figure reveals the nuanced variations in temperature profiles over time, demonstrating 181 the explicit scheme’s capability to provide a detailed representation of the system’s 182 behavior. 183 Contrastingly, Figure 3 presents the temporal evolution of temperature distribution 184 using an implicit scheme with different time steps. Each subplot offers a comprehensive 185 comparison of the model’s performance under distinct temporal resolutions, revealing 186 the implicit scheme’s versatility. Even in the case of larger time steps, such as ∆ t= 10 187 s (Figure 3d), the implicit scheme demonstrates convergence, underscoring its stability 188 across a broader temporal range. 189 6 The comparison between Figure 2 and Figure 3 emphasizes that the implicit scheme 190 maintains its superiority even when considering single steps over the specified time 191 range. The implicit scheme’s ability to handle larger time steps without compromising 192 convergence or accuracy positions it as a robust and efficient choice for simulating the 193 1D Heat Equation, particularly in scenarios where computational resources or time 194 constraints necessitate the use of larger time steps. 195 Fig. 2 : Temporal evolution of temperature distribution utilizing explicit scheme with time step ∆ t= 0.00025 s. The figure illustrates the dynamic changes in temperature over time, providing a detailed insight into the computational model’s behavior under the specified temporal resolution. 7 (a) (b) (c) (d) Fig. 3 : Temporal evolution of temperature distribution utilizing the implicit Scheme with different time steps. (a) ∆ t= 0.1 s, (b) ∆ t= 0.25 s, (c) ∆ t= 2.5 s, and (d) ∆ t= 10 s. These figures elucidate the dynamic changes in temperature over time, providing comprehensive insights into the computational model’s behavior under distinct temporal resolutions for the 1D Heat Equation. In evaluating the performance of our 1D Heat Equation simulation employing both 196 explicit and implicit methods, we undertook a comparative analysis of runtime met- 197 rics with and without the integration of the Numba library. The explicit method, 198 when executed without Numba, exhibited an average total runtime of 8.324 s, with a 199 standard deviation of 0.593 seconds and a range spanning from 6.681 to 9.602 s. In 200 8 contrast, the explicit method leveraging Numba demonstrated a substantial enhance- 201 ment, showcasing an average runtime of 4.035 s, a diminished standard deviation of 202 0.251 s, and a narrower range between 3.524 and 4.744 s. This observation suggests a 203 notable acceleration in the overall computational efficiency facilitated by Numba. 204 Turning our attention to the implicit method, its execution without Numba 205 resulted in an average total runtime of 9.970 s, accompanied by a standard deviation 206 of 0.967 s and a range varying from 7.291 to 12.054 s. The incorporation of Numba into 207 the implicit method led to a remarkable improvement, yielding an average runtime of 208 1.195 s, a reduced standard deviation of 0.215 s, and a narrower range spanning from 209 0.758 to 1.802 s. The most noteworthy acceleration in runtime was observed in the 210 implicit method, indicating that Numba has a particularly pronounced impact on this 211 scheme. 212 The effectiveness of Numba can be attributed to its ability to translate Python code 213 directly into machine language via LLVM, thereby bypassing the interpretative over- 214 head associated with Python. This direct translation contributes significantly to the 215 efficiency gains observed in both explicit and implicit methods. Furthermore, the paral- 216 lelization process occurring behind the scenes, facilitated by Numba, further enhances 217 the computational speed. In the context of the implicit method, which involves the 218 recursive solution of a tridiagonal matrix structure, the inherent parallelizability of 219 such tasks results in a more substantial reduction in total runtime. 220 The visual representation of these findings is presented in Figure 4 and Figure 5, 221 illustrating the stark differences in total runtime between control experiments (with- 222 out Numba) and test experiments (with Numba) for the explicit and implicit schemes, 223 respectively. The visualizations align with the quantitative results, providing a com- 224 prehensive illustration of the efficiency gains achieved through the integration of 225 Numba. 226 These descriptive statistical results demonstrate that Numba significantly accel- 227 erates the total runtime of 1D Heat Equation simulations across both explicit and 228 implicit methods. While notable improvements were observed in both schemes, the 229 most rapid enhancement occurred in the implicit method. This highlights the potential 230 of Numba as a valuable tool for optimizing Python-based numerical computing codes, 231 particularly in scenarios involving intricate computations and recursive structures. 232 9 (a) (b) (c) (d) Fig. 4 : Total runtime for the 1D Heat Equation’s explicit scheme with Numba (test, red) and the standard implementation (control, blue) over 100 runs, showcasing (a) distribution plots, (b) boxplots, (c) violin plots, and (d) ECDFs, rendered using Mat- plotlib and seaborn [30]. 10 (a) (b) (c) (d) Fig. 5 : Same as Figure 4, but for the implicit scheme. The examination of results is underpinned by a meticulous analysis of discernible 233 differences between the control and test experiments. The vivid portrayal of these 234 distinctions through distribution boxplots, violin plots, and Empirical Cumulative 235 Distribution Functions (ECDFs) in Figures 4 and 5 offers a rich visual narrative of 236 the explicit and implicit schemes’ performance variations. 237 To rigorously evaluate the statistical significance of these disparities, we employed 238 robust statistical tests—specifically, the Mann-Whitney U test and the Wilcoxon 239 Signed Rank test. The Mann-Whitney U test yielded a consistent U statistics score 240 of 10,000 for both explicit and implicit schemes, accompanied by an impressively low 241 11 p-value of <0.01. Simultaneously, the Wilcoxon Signed Rank test consistently pro- 242 duced a W-statistics score of 0, again with a p-value <0.01. These findings not only 243 indicate statistical significance but also emphasize the reliability and stability of the 244 observed differences. 245 The persistent U statistics scores of 10,000 for both schemes, coupled with identical 246 p-values, underscore the robustness of our results across different numerical meth- 247 ods. This uniformity in statistical measures fortifies the argument for the consistent 248 superiority of the optimized schemes. 249 The W-statistics score of 0 further substantiates the claim of statistical significance, 250 emphasizing the reliability and stability of results across both explicit and implicit 251 schemes. This score signifies a high level of consistency in performance improvements, 252 adding a layer of confidence to the conclusion that the observed enhancements are not 253 sporadic but systematically linked to the incorporation of Numba. 254 In light of these statistically robust measures, we assert with confidence that the 255 incorporation of Numba yields a statistically significant optimization in accelerating 256 the numerical solution of the 1D Heat Equation. The p-values below 0.01 underline the 257 robustness of the results, reinforcing the assertion that the observed improvements are 258 not coincidental but are indeed attributable to the deliberate optimization strategy. 259 This high level of statistical significance enhances the credibility of our findings and 260 underscores the effectiveness of leveraging Numba for augmenting the computational 261 efficiency of numerical methods applied to the 1D heat equation. 262 4 Conclusion 263 In conclusion, our exploration into optimizing numerical solvers for the 1D Heat 264 Equation using Numba has yielded significant advancements in computational effi- 265 ciency. The interplay between Python’s elegance and Numba’s JIT compilation 266 capabilities has allowed us to craft a bespoke solver that not only encapsulates the 267 readability and versatility of Python but also harnesses the computational alacrity 268 afforded by Numba. Through a granular presentation of our Numba-based solver, 269 we navigated the complexities of implementation, optimization strategies, and rigor- 270 ous performance evaluations. The explicit and implicit finite difference methods were 271 examined, and the efficiency gains achieved through Numba integration were quan- 272 tified and statistically validated. Our results demonstrate that Numba substantially 273 accelerates the total runtime of 1D Heat Equation simulations. The explicit method 274 showcased an average runtime reduction from 8.324 s to 4.035 s, while the implicit 275 method exhibited an even more pronounced improvement, decreasing from 9.970 s to 276 1.195 s. Statistical tests, including the Mann-Whitney U test and the Wilcoxon Signed 277 Rank test, consistently confirmed the statistical significance of these enhancements. 278 Moving forward, the success of our Numba-based solver opens avenues for future 279 research and development in several directions. Firstly, we aim to extend the optimiza- 280 tion strategies to handle multidimensional heat equations, investigating the challenges 281 and opportunities presented by higher-dimensional systems and exploring how Numba 282 can be leveraged for efficient solutions. Additionally, we plan to further explore and 283 12 optimize parallelization techniques within Numba to exploit parallel hardware archi- 284 tectures fully. This involves investigating the impact of parallelization on both explicit 285 and implicit methods, considering the scalability to larger problem sizes. 286 Furthermore, we intend to implement dynamic parameter optimization strategies 287 within Numba to adaptively adjust parameters such as time step sizes for different 288 simulation scenarios. This could enhance the solver’s versatility and applicability to 289 a broader range of problems. Moreover, we envision the integration of Numba-based 290 solvers with advanced numerical techniques or machine learning methods for enhanced 291 accuracy and efficiency. Exploring hybrid approaches that combine the strengths of 292 different solvers will be a key focus. 293 Finally, the application of the optimized solver to real-world problems in diverse 294 fields such as materials science, electronics, and environmental modeling is a critical 295 next step. Collaborating with domain experts to validate the solver’s performance in 296 practical settings will provide valuable insights and contribute to the solver’s real-world 297 impact. 298 In summary, our Numba-based solver lays the groundwork for a new era in efficient 299 numerical solutions for the 1D Heat Equation. As we delve deeper into multidimen- 300 sional problems, explore advanced parallelization, and apply our solver to real-world 301 challenges, the synergy between Python and Numba promises to reshape the landscape 302 of scientific computing. 303 Acknowledgments. This work was supported by the Dean’s Distinguished Fel- 304 lowship at the University of California, Riverside (UCR) 2023 and ITB Research, 305 Community Service and Innovation Program (PPMI-ITB) 2023. The code and total 306 runtime dataset for this paper are hosted on our GitHub: https://github.com/ 307 sandyherho/1DHeatEqnNumba. 308