# A Truly Second-Order and Unconditionally Stable Thermal Lattice Boltzmann Method
## 4 Numerical Tests
### 2D Natural Convection in a Square Cavity  
Natural convection in a square cavity is a classic benchmark test to validate numerical methods for thermal flows. In this section, we use the example in two-dimensional scenario to assess accuracy, convergence, and stability performance of USTLBM. Physical configuration of this problem is shown in Figure 2. Key non-dimensional parameters that affect the flow pattern are the Rayleigh number and the Prandtl number defined as  
\[Ra=\frac{g\beta\cdot\Delta T\cdot L^{3}}{v\cdot\chi}=\frac{V_{c}^{2}\cdot L^{ 2}}{v\cdot\chi} \tag{66}\]  
\[Pr=\frac{v}{\chi} \tag{67}\]  
where \(\Delta T\) is the temperature difference between the left hot wall and the right cold wall; and \(V_{c}\) is the characteristic velocity defined as \(V_{c}=\sqrt{g\beta L\cdot\Delta T}\). In our simulations, \(Pr\) is set to be 0.71, and \(V_{c}=1\) to fulfill the low Mach number assumption. Various Nusselt numbers are defined to quantify the heat transfer rate. Specifically, the volume averaged Nusselt number and the average Nusselt number on the hot wall are  
\[\overline{Nu}=\frac{L}{\chi\cdot\Delta T}\cdot\frac{1}{L^{2}}\underset{\Omega }{\iint}\left(uT-\chi\frac{\partial T}{\partial x}\right)\text{d}\Omega \tag{68}\]  
\[Nu_{0}=\frac{1}{L\cdot\Delta T}\underset{x=0}{\int}\frac{\partial T}{\partial x }\text{d}y \tag{69}\]  
Initially, we use the case with \(Ra=10^{4}\) to validate spatial accuracy of USTLBM. Since analytical solution to this problem is unavailable, the benchmark result of the volume averaged Nusselt number obtained by the mesh size of 301 \(\times\) 301 is used as the reference data \(\overline{Nu}_{a}\). Various mesh sizes of 11 \(\times\) 11, 21 \(\times\) 21, 41 \(\times\) 41, 61 \(\times\) 61 and 81 \(\times\) 81 are adopted and corresponding relative errors are calculated through the following formula  
\[E=\frac{\left|\overline{Nu}-\overline{Nu}_{a}\right|}{\left|\overline{Nu}_{a }\right|} \tag{70}\]  
The relative errors obtained at different mesh sizes are depicted in Figure 3. As can be seen, the linearly fitted line of the data reveals a slope of 1.956 in log scale, which indicates that USTLBM can ensure the second-order of accuracy in space. The result obtained by STLBM is also included in this figure. The fitted line of the data given by STLBM also has a slope that is very near to 2.0, which validates that STLBM can ensure the second-order of accuracy in space. However, it should be noted that the order of accuracy and the numerical error (e.g., the relative error to the reference data in this example) are two different concepts. The numerical error is usually written as \(\text{Error}=O(\mu^{n})=C\cdot h^{n}\), where \(h\) is the mesh spacing, \(n\) is the order of accuracy, and \(C\) is a constant. Clearly, the order of accuracy is the slope of straight line of log (error) vs. \(\log(h)\). Even if the order of accuracy \(n\) is the same, the relative numerical error of one method can still be smaller than that of another method at the same mesh size (i.e., \(C\) value of two methods may not be the same). Detailed comparisons show that, at the same mesh spacing, relative error of the results given by USTLBM is lower than that obtained by STLBM. Therefore, it is reasonable to expect more accurate results at the same mesh resolution when applying USTLBM. More cases in high Rayleigh numbers of \(Ra=10^{5}\), \(10^{6}\), \(10^{7}\) and \(10^{8}\) are simulated by both USTLBM and STLBM with different mesh resolutions. Note that, in cases of \(Ra=10^{7}\) and \(10^{8}\), non-uniform meshes are adopted to capture severe changes of flow fields near boundaries. Lagrange interpolation algorithm is applied to determine functional values on locations where the streaming process initiates. Representative flow parameters are tabulated in Tables 1-4. It is noteworthy that, in all cases, both USTLBM and STLBM do not diverge at very coarse mesh size of \(11\times 11\), and reasonable results of flow parameters can be achieved.