Lattice Boltzmann Simulation of Coupling Heat Transfer between Solid and Gas Phases of Nanoporous Materials

In order to deeply study the heat conduction of nanoporous aerogel, a model of gas-solid heat conduction was established based on the microstructure of aerogel. The model was divided into two subdomains with uniform mesh because of the different gas-solid characteristics, and simulation was performed on each domain using the lattice Boltzmann method. The value of temperature on the boundaries of subdomains was determined by interpolation. Finally, the temperature distribution and the thermal conductivity were maintained. It can be concluded that when the gas-phase scale was fixed, the temperature distribution of the solid phase became more uniform when the scale increased; when the solid-phase scale was fixed, the temperature jump on the gas-solid interface decreased with the increase in the gas-phase scale; and the thermal conductivity of gas-solid coupling varied with the scale of the gas phase or solid phase, showing a scale effect in varying degrees.


Introduction
Aerogels, as a super insulation material, have a unique structure that is an open nanoporous and continuous three-dimensional skeleton. Therefore, they exhibit outstanding properties such as high porosity, high specific surface area and low thermal conductivity [1]. They have been applied in many fields such as aerospace, industrial insulation and military technology [2]. In the thermal transport of aerogel, the heat needs to continuously pass through the porous network composed of gas and solid. The complex microscopic characteristics and the scale effect of gas and solid make it difficult to study the thermal conductivity. With the depletion of natural resources, all countries in the world are committed to the conservation and efficient use of energy. The insulation requirements of aerogel insulation materials are becoming more and more urgent. Accordingly, it is significant to further study the gas-solid coupling thermal conductivity of aerogels in order to optimize and improve their adiabatic performance.
At present, the researches on the gas-solid coupling heat transfer of aerogel mainly focus on theoretical and numerical calculation methods. In terms of theoretical studies, Swimm et al. proposed a calculation model of coupling the thermal conductivity of aerogels and studied the coupling effect of their effective thermal conductivity [3]. Based on the microstructure of aerogels and the model, Zhao J et al. improved the prediction accuracy of the thermal conductivity of aerogels [4]. Given that the above two models are based on the thermal bridge effect, the calculation process is extremely difficult. In addition, Wei et al. used a cubic array composed of small spheres to establish a gas-solid coupling heat conduction model for aerogels [5]. The calculation results showed that the nanostructure and the nanoscale effect of the solid grains were the main reasons for the extremely low thermal conductivity of the materials. In terms of numerical methods, Coquard R. et al. established a model for the overall heat transfer inside aerogels with cellulose [6]. Under the assumption of "uniform approximation", the heat-transfer law of gas-solid two phases was studied. Zhao J. et al. established a three-dimensional numerical model of silica aerogel based on the random aggregation structure of porous secondary nanoparticles and studied the overall thermal conductivity of aerogel by considering the coupling thermal conductivity [7].
In recent years, the lattice Boltzmann method (LBM) has been widely used in the field of micro-nanoscale heat transfer [8]. It can not only reveal the heat transfer and the intrinsic mechanism of porous insulating materials, but it can also make up for the shortcomings of the macro methods. In addition, the LBM is high efficiency and stable. It has become a powerful tool for the analysis of micro-scale heat transfer. Song proposed a local-effective-viscosity multirelaxation-time lattice Boltzmann pore-network coupling model (LEV-LBM-PNM) to simulate gas transport in nanoporous media [9]. The proposed LEV-LBM-PNM accurately predicts gas-apparent permeability by accounting for gas slip in irregular pore shapes and surface roughness. Qu established a modified lattice LBM to predict effective thermal conductivity of aerogel materials [10]. The modified LBM scheme introduces an additional coefficient to regulate significant differences in inherent thermal conductivity between solid and gas phases. Kan proposed a random internal morphology and structure generation-growth method [11], The model was imported into the lattice Boltzmann algorithm as a fully resolved geometry and used to investigate the effects on heat transfer in nanoporous materials.
In this paper, the gas-solid coupling heat conduction model is established, and the temperature distribution and thermal conductivity of silica-aerogel gas-solid coupling are simulated by LBM based on the region splitting. The thermal transfer of gas-solid coupling heat conduction is analyzed and discussed in depth. Figure 1 displays the model of gas-solid coupling heat conduction of aerogels. Figure 1a,b are a three-dimensional model and a two-dimensional gas-solid center section, respectively. The model consists of gas-solid two phases, which are air and silica, respectively, and the characteristic scales along the x direction are expressed by d g and d s . The length, width and height of the model were expressed as L x , L y , L z . The red cross section is the interface of the gas-solid two phases. The front and back walls are kept at a constant temperature; the temperature of the front wall is T 1 = 300 K, and the temperature of the back wall is T 2 = 301 K; the other walls were adiabatic. The heat flow transferred from the high temperature wall (front wall) to the low temperature wall (back wall) along the x-axis. Nanomaterials 2022, 12, 3424 2 of cellulose [6]. Under the assumption of "uniform approximation", the heat-transfer law gas-solid two phases was studied. Zhao J. et al. established a three-dimensional numeri model of silica aerogel based on the random aggregation structure of porous seconda nanoparticles and studied the overall thermal conductivity of aerogel by considering t coupling thermal conductivity [7]. In recent years, the lattice Boltzmann method (LBM) has been widely used in the fie of micro-nanoscale heat transfer [8]. It can not only reveal the heat transfer and the intr sic mechanism of porous insulating materials, but it can also make up for the shortco ings of the macro methods. In addition, the LBM is high efficiency and stable. It has b come a powerful tool for the analysis of micro-scale heat transfer. Song proposed a loc effective-viscosity multirelaxation-time lattice Boltzmann pore-network coupling mod (LEV-LBM-PNM) to simulate gas transport in nanoporous media [9]. The proposed LE LBM-PNM accurately predicts gas-apparent permeability by accounting for gas slip irregular pore shapes and surface roughness. Qu established a modified lattice LBM predict effective thermal conductivity of aerogel materials [10]. The modified LBM schem introduces an additional coefficient to regulate significant differences in inherent therm conductivity between solid and gas phases. Kan proposed a random internal morpholo and structure generation-growth method [11], The model was imported into the latt Boltzmann algorithm as a fully resolved geometry and used to investigate the effects heat transfer in nanoporous materials.

Physical Model
In this paper, the gas-solid coupling heat conduction model is established, and t temperature distribution and thermal conductivity of silica-aerogel gas-solid coupling a simulated by LBM based on the region splitting. The thermal transfer of gas-solid co pling heat conduction is analyzed and discussed in depth. Figure 1 displays the model of gas-solid coupling heat conduction of aerogels. Figu 1a,b are a three-dimensional model and a two-dimensional gas-solid center section, spectively. The model consists of gas-solid two phases, which are air and silica, resp tively, and the characteristic scales along the x direction are expressed by dg and ds. T length, width and height of the model were expressed as Lx, Ly, Lz. The red cross section the interface of the gas-solid two phases. The front and back walls are kept at a consta temperature; the temperature of the front wall is T1 = 300 K, and the temperature of t back wall is T2 = 301 K; the other walls were adiabatic. The heat flow transferred from t high temperature wall (front wall) to the low temperature wall (back wall) along the axis.

Gas-Phase Basic Algorithm
The three-dimensional fifteen-speed (D3Q15) lattice Boltzmann model was used to solve the gas steady heat conduction of aerogel. The evolution equation was as follows: where subscripts g, k were the gas phase and the discrete velocity directions of the lattice point, respectively; g r was the space position vector , k c was the discrete lattice velocity (as shown in Figure 2); g t was the time, where T was temperature; ρ was the air density; cp was the specific heat capacity of the air.
where λ was the thermal conductivity of the air; c was the lattice speed that theoretically could take any positive value only to insure the τ value within (0.5,2); c = δx/δt was the The three-dimensional fifteen-speed (D3Q15) lattice Boltzmann model was used to solve the gas steady heat conduction of aerogel. The evolution equation was as follows: where subscripts g, k were the gas phase and the discrete velocity directions of the lattice point, respectively; r g was the space position vector, c k was the discrete lattice velocity (as shown in Figure 2); t g was the time, ∆t g was the time step; τ was the dimensionless relaxation time; f k with f eq k were, respectively, the internal energy distribution function and the corresponding equilibrium distribution function of the discrete velocity directions of the different lattice points in the gas phase.
where T was temperature; ρ was the air density; c p was the specific heat capacity of the air.
where λ was the thermal conductivity of the air; c was the lattice speed that theoretically could take any positive value only to insure the τ value within (0.5,2); c = δx/δt was the space step [12,13]. Temperature T and heat flux q could be obtained, respectively, according to the following formula [14]: When the temperature reached equilibrium state, the effective thermal conductivity could be calculated on the temperature T and heat flux q.

Solid-Phase Basic Algorithm
The three-dimensional fifteen-speed (D3Q15) lattice Boltzmann model was used to solve the solid steady heat conduction of aerogel. The evolution equation was as follow: where subscripts s, k were the solid phase and the discrete velocity directions of the lattice point, respectively; s r was the space position vector , k v was the discrete lattice velocity (as shown in Figure 2); s t was the time, s t Δ was the time step; β was a weight function; k e with eq k e were, respectively, the phonon energy density distribution function and the equilibrium phonon energy density distribution function of the discrete velocity directions of the different lattice points in the solid phase. The total phonon energy density e is calculated in reference [15]. The total energy density e of each lattice point is discretized into the energy density k e of each discrete direction.
Since the phonon energy density e and the equilibrium phonon energy density eq k e of all lattice points need to be redefined at each time step, assuming that the probability of phonon scattering in all directions is equal, i.e. isotropic, the equilibrium phonon energy density When the temperature reached equilibrium state, the effective thermal conductivity could be calculated on the temperature T and heat flux q.

Solid-Phase Basic Algorithm
The three-dimensional fifteen-speed (D3Q15) lattice Boltzmann model was used to solve the solid steady heat conduction of aerogel. The evolution equation was as follow: where subscripts s, k were the solid phase and the discrete velocity directions of the lattice point, respectively; r s was the space position vector, v k was the discrete lattice velocity (as shown in Figure 2); t s was the time, ∆t s was the time step; β was a weight function; e k with e eq k were, respectively, the phonon energy density distribution function and the equilibrium phonon energy density distribution function of the discrete velocity directions of the different lattice points in the solid phase.
The total phonon energy density e is calculated in reference [15]. The total energy density e of each lattice point is discretized into the energy density e k of each discrete direction.
Since the phonon energy density e and the equilibrium phonon energy density e eq k of all lattice points need to be redefined at each time step, assuming that the probability of phonon scattering in all directions is equal, i.e. isotropic, the equilibrium phonon energy density where the lattice speed v = ∆r s /∆t s , ∆r s was the solid-phase space step. Weight function β = ∆t s /τ 0 , τ 0 was the particle relaxation time.
The normal heat-transfer rate of the yoz plane within the model is Nanomaterials 2022, 12, 3424 5 of 13

Boundary Conditions
In the LBM, for calculation accuracy, stability and computational efficiency, the boundary conditions should be reasonably selected according to the actual situation. Since the influences of solid skeleton were not taken into account in the gas heat conduction, the mirror bounce format was adopted for no friction losses. Moreover, the energy distribution function could free transfer in this boundary condition. That was equivalent to removing the boundary, and the results would not be affected by the width L y and the height L z . The L x direction is the constant temperature boundaries T 1 and T 2 . The specific forms were implemented as follow: where r b was the boundary lattice point, r f = r b − e i δt were internal lattice point; f i was a mirror-symmetric distribution function of f i . The f i could be obtained by the inner lattice points that were adjacent the boundary wall.

Calculation Method of Effective Thermal Conductivity
Owing to both the gas-phase region and the solid-phase region being stable and the heat flow in the two sub-regions being equal, q is used to represent the heat flow in the whole region, that is q = q s = q g . Because the characteristic scales of the whole system region are smaller than or close to the average free path of the gas molecules and phonons, obvious scale effects appear in both sub-regions, and the temperature distribution of materials is no longer continuous. The effective thermal conductivity can be obtained according to Fourier's law and the geometric relationship in Figure 3.
where q is the system heat flow; ∆T a is the jumping for the gas-phase boundary temperature; ∆T c is the jumping for the solid-phase boundary temperature.

Region Splitting
Although the hybrid LBM has been proven to provide great advantages over the structured LB in terms of geometry flexibility and the treatment of interfaces [16,17], the focus of this paper is to study the performance of gas-solid coupling heat transfer at the

Region Splitting
Although the hybrid LBM has been proven to provide great advantages over the structured LB in terms of geometry flexibility and the treatment of interfaces [16,17], the focus of this paper is to study the performance of gas-solid coupling heat transfer at the interface. Hence, considering the accuracy and simplicity, the "region splitting" is selected.
The basic idea of region splitting is to divide the whole system region into different sub-regions, and then using the corresponding LBM model to calculate each sub-region independently [18].
According to the gas-solid coupled thermal conduction model, the gas-solid coupled thermal conduction model is divided into two sub-regions: the gas phase and the solid phase, which are represented by Ω g and Ω s , Γ as = Ω a ∩ Ω s , because the average free path of a gas molecule (about 70 nm) is much larger than that of a phonon (about 0.6 nm) of silica solid. To simulate the heat transfer at the solid boundary, we set a virtual boundary node at the gas region XA, as shown in Figure 4. Therefore, different size meshes are used to subdivide the sub-regions, which are set up as ∆x a = ∆y a = ∆z a , ∆x s = ∆y s = ∆z s . Generally, the time step and space step in the gas-phase and solid-phase regions are different, so the physical quantities involved in the virtual boundary need to be calculated by special methods.  Set / g s n x x = Δ Δ as an integer, that is, the lattice size of the gas phase is several times larger than that of solid phase. Because the lattice size of gas phase and solid phase is different, some lattices will coincide on the virtual boundary, and some lattices will disperse on the virtual boundary, as shown in Figure 5. Because the lattice size of the gassolid two-phase is different, the time step of the gas-solid two-phase also has a special relationship. If the gas-phase region evolves from time 0 t to time The physical quantities on lattice A at the virtual boundary time t are calculated by the second-order interpolation method to ensure the second-order accuracy of the calculation results. Since there is no temperature gradient in the y and z directions, and the boundary is slip free adiabatic, the specific implementation steps are as follows: (1) If g s A∈Ω Ω  , as shown in Figure 5a, then ( ) ( ) ( )  Set n = ∆x g /∆x s as an integer, that is, the lattice size of the gas phase is several times larger than that of solid phase. Because the lattice size of gas phase and solid phase is different, some lattices will coincide on the virtual boundary, and some lattices will disperse on the virtual boundary, as shown in Figure 5. Because the lattice size of the gas-solid twophase is different, the time step of the gas-solid two-phase also has a special relationship. If the gas-phase region evolves from time t 0 to time t g = t 0 + ∆t g , the solid-phase region needs to evolve n time steps: t s = t 0 + ∆t s , t s = t 0 + 2∆t s , · · · , t s = t 0 + n∆t s . The physical quantities on lattice A at the virtual boundary time t are calculated by the second-order interpolation method to ensure the second-order accuracy of the calculation results. Since there is no temperature gradient in the y and z directions, and the boundary is slip free adiabatic, the specific implementation steps are as follows: (1) If A ∈ Ω g ∩ Ω s , as shown in Figure 5a, then Here C ∈ Ω g and B ∈ Ω s are the nearest lattice points in the vertical direction of the virtual boundary.
(2) If A / ∈ Ω g ,A ∈ Ω s , as shown in Figure 5b, then Here D ∈ Ω g , α = |XC|/|DC|.  The gas-solid two-phase coupling is performed by the second-order interpolation method. According to the region splitting, the gas-solid two-phase sub-regions are calculated by their respective LBM models, and the coupling calculation of the two regions is made by counting the macroscopic physical parameters, T, at the interface. Finally, the heat transport process of the whole gas-solid is obtained. Figure 5 shows the numerical simulation results of the thermal conductivity of the gas-phase region at the pore scale of 10-80 nm, which is compared with the Zeng model [19] under the same conditions. It can be seen from the figure that the simulation results of this paper are basically consistent with the variation of the theoretical analysis results of the Zeng model, which is slightly smaller than the Zeng model. This is because the potential was weak, and it only considered pairwise interaction between atoms. As a result, the simulation results deviated from this paper in 60~70 nm, although they have the same change trend. Hence, the reliability of the LBM program could be proved. Figure 6 shows a comparison between the simulation results of thermal conductivity of the solid-phase region and the results of reference [15] under the same conditions. It can be seen from the figure that the simulation results of the solid-phase region in this paper are basically consistent with the variation trend of reference [15]. The gas-solid two-phase coupling is performed by the second-order interpolation method. According to the region splitting, the gas-solid two-phase sub-regions are calculated by their respective LBM models, and the coupling calculation of the two regions is made by counting the macroscopic physical parameters, T, at the interface. Finally, the heat transport process of the whole gas-solid is obtained. Figure 5 shows the numerical simulation results of the thermal conductivity of the gas-phase region at the pore scale of 10-80 nm, which is compared with the Zeng model [19] under the same conditions. It can be seen from the figure that the simulation results of this paper are basically consistent with the variation of the theoretical analysis results of the Zeng model, which is slightly smaller than the Zeng model. This is because the potential was weak, and it only considered pairwise interaction between atoms. As a result, the simulation results deviated from this paper in 60~70 nm, although they have the same change trend. Hence, the reliability of the LBM program could be proved. Figure 6 shows a comparison between the simulation results of thermal conductivity of the solid-phase region and the results of reference [15] under the same conditions. It can be seen from the figure that the simulation results of the solid-phase region in this paper are basically consistent with the variation trend of reference [15]. When the characteristic scale of the solid-phase region is between 0 and 10 nm, the thermal conductivity shows a significant jump, showing the scale-effect characteristics, while the thermal conductivity is close to 1.47 W·m −1 ·K −1 of silica in the macroscopic state at 10-60 nm. Figure 7 is a steady-state temperature distribution cloud diagram of a two-dimensional gas-solid center section, in which the gas-phase region d g is 60 nm, and the solid-phase region d s is 2 nm, 3 nm, 4 nm, and 5 nm, respectively. Because the solid-phase scale is much smaller than the gas-phase scale, the solid-phase region is presented in the form of enlarged views. From the temperature distribution profiles in Figure 7a-d, since the gas-phase scale is 60 nm, it can be seen that the gas-phase characteristic scale is smaller than the average free path of the gas molecule and the discontinuous temperature diffusion mainly concentrated at the boundary of x = 0 nm and x = 60 nm, while the temperature distribution in the central region is relatively uniform. From the corresponding temperature distribution diagram of the solid-phase region, it can be seen that the scattering between phonons is dominant because the characteristic scale of the solid phase is much larger than the average free path of the silica phonons at 0.6 nm. At the same time, the characteristics of the concentrated distribution of the solid-phase boundary temperature are weakened. With the increase in the characteristic scale of the solid phase, the transition between the different parts of the temperature distribution becomes more uniform and tends to the temperature diffusion phenomenon under the macroscopic state.    Figure 8 shows the two-dimensional steady-state temperature distribution of the g solid central section of the gas-solid coupled heat conduction model. It can be seen fr the graph that the gas-phase region dg is 60 nm, and since the gas-phase characteristic sc is small, the free movement of the gas molecules is severely restricted, which hinders thermal transport of the heat molecules to the cold molecules, and a strong boundary s tering effect occurs. Therefore, the temperature jump is obvious at the boundary x = 0 and x = 60 nm. It can be seen from the large graph of the respective solid-phase region t there is a significant temperature jump at the interface of gas-solid coupling x = 60 n and The temperature jump of the gas-phase region is significantly higher than the t perature jump of the solid-phase region. This is because the characteristic scale of solid-phase region is much larger than the average free path of the silica phonons which time the thermal transport between the phonons plays a dominant role, and scattering between the phonons and the boundaries is relatively weakened. Hence, temperature jump between the two boundaries is weakening.  Figure 8 shows the two-dimensional steady-state temperature distribution of the gassolid central section of the gas-solid coupled heat conduction model. It can be seen from the graph that the gas-phase region d g is 60 nm, and since the gas-phase characteristic scale is small, the free movement of the gas molecules is severely restricted, which hinders the thermal transport of the heat molecules to the cold molecules, and a strong boundary scattering effect occurs. Therefore, the temperature jump is obvious at the boundary x = 0 nm and x = 60 nm. It can be seen from the large graph of the respective solid-phase region that there is a significant temperature jump at the interface of gas-solid coupling x = 60 nm, and The temperature jump of the gas-phase region is significantly higher than the temperature jump of the solid-phase region. This is because the characteristic scale of the solid-phase region is much larger than the average free path of the silica phonons, at which time the thermal transport between the phonons plays a dominant role, and the scattering between the phonons and the boundaries is relatively weakened. Hence, the temperature jump between the two boundaries is weakening.   Figure 9 is a temperature jump diagram at the gas-solid interface when the solidphase region is 5 nm. It can be seen from the graph that the temperature jump at the interface of the gas-solid coupling decreases with the increase in the gas-phase scale. The temperature jump at 10-70 nm decreases sharply when it reaches 70 nm, and the decrease amplitude begins to be gentle.    Figure 9 is a temperature jump diagram at the gas-solid interface when the solid-phase region is 5 nm. It can be seen from the graph that the temperature jump at the interface of the gas-solid coupling decreases with the increase in the gas-phase scale. The temperature jump at 10-70 nm decreases sharply when it reaches 70 nm, and the decrease amplitude begins to be gentle.   Figure 9 is a temperature jump diagram at the gas-solid interface when the solidphase region is 5 nm. It can be seen from the graph that the temperature jump at the interface of the gas-solid coupling decreases with the increase in the gas-phase scale. The temperature jump at 10-70 nm decreases sharply when it reaches 70 nm, and the decrease amplitude begins to be gentle.  This is because when the gas-phase scale is smaller than the mean free path of the gas molecules, the boundary scattering effect is obvious, and the temperature jump at the boundary will be more violent. With the increase in the gas-phase scale, the discontinuous temperature distribution at the boundary caused by the boundary scattering effect will gradually weaken at same time, leading to the decrease in the boundary temperature jump value. Figure 10 is a schematic diagram of the dimensionless temperature distribution of the gas-solid two-phase system. As shown in Figure 10a, since the gas-phase characteristic scale is 60 nm, the temperature distributions in the dimensionless gas-phase regions in the four cases are in good agreement, and the boundary temperature jumps are consistent; as shown in Figure 10b, the temperature jumps in the solid-phase region decrease with the increase in the solid-phase scale, and the temperature diffusion is more and more close to the macroscopic heat transport characteristics. This is because the solid-phase scale is larger than the average free path of phonon 0.6 nm. With the further increase in the solid-phase scale, the scattering effect between the phonon and the boundary will gradually disappear, the scattering between phonon and phonon will be enhanced, and the energy transfer will be smoother. Hence, the macroscopic thermal transport characteristics will gradually appear. This is because when the gas-phase scale is smaller than the mean free path of the ga molecules, the boundary scattering effect is obvious, and the temperature jump at th boundary will be more violent. With the increase in the gas-phase scale, the discontinuou temperature distribution at the boundary caused by the boundary scattering effect wi gradually weaken at same time, leading to the decrease in the boundary temperature jum value. Figure 10 is a schematic diagram of the dimensionless temperature distribution o the gas-solid two-phase system. As shown in Figure 10a, since the gas-phase characterist scale is 60 nm, the temperature distributions in the dimensionless gas-phase regions in th four cases are in good agreement, and the boundary temperature jumps are consistent; a shown in Figure 10b, the temperature jumps in the solid-phase region decrease with th increase in the solid-phase scale, and the temperature diffusion is more and more close t the macroscopic heat transport characteristics. This is because the solid-phase scale larger than the average free path of phonon 0.6 nm. With the further increase in the solid phase scale, the scattering effect between the phonon and the boundary will graduall disappear, the scattering between phonon and phonon will be enhanced, and the energ transfer will be smoother. Hence, the macroscopic thermal transport characteristics wi gradually appear.

Simulation Results and Analysis
(a) Gaseous region (b) Solid region Figure 10. Dimensionless gas-solid temperature distribution. * dimensionless parameters. Figure 11 shows the relationship between the gas-solid coupling thermal conductiv ity and the solid-gas scale. It can be seen from Figure 11a that when the solid-phase sca is less than 1 nm, the increase in thermal conductivity is small, and when the solid-phas scale is greater than 1 nm, the thermal conductivity begins to increase rapidly. As can b seen from Figure 11b, when the gas-phase scale is between 10 and 70 nm, the increase i thermal conductivity is larger, and when the gas-phase scale reaches about 70 nm, th change trend of thermal conductivity begins to be stable. This is because the characteristic scales of the gas-solid two phases are smaller o closer to the mean free path of their molecules and phonons, and boundary scattering wi severely restrict the heat transport between the molecules and phonons, showing the ob vious scale effect. With the increase in scale, scattering will rapidly weaken and therma conductivity will correspondingly increase.  Figure 11 shows the relationship between the gas-solid coupling thermal conductivity and the solid-gas scale. It can be seen from Figure 11a that when the solid-phase scale is less than 1 nm, the increase in thermal conductivity is small, and when the solid-phase scale is greater than 1 nm, the thermal conductivity begins to increase rapidly. As can be seen from Figure 11b, when the gas-phase scale is between 10 and 70 nm, the increase in thermal conductivity is larger, and when the gas-phase scale reaches about 70 nm, the change trend of thermal conductivity begins to be stable. This is because the characteristic scales of the gas-solid two phases are smaller or closer to the mean free path of their molecules and phonons, and boundary scattering will severely restrict the heat transport between the molecules and phonons, showing the obvious scale effect. With the increase in scale, scattering will rapidly weaken and thermal conductivity will correspondingly increase.

Conclusions
In this paper, the gas-solid coupling heat conduction model of nanoporous aerogels is established. Based on the theory of regional splitting, the gas-solid two-phase simulation is carried out by using the respective LBM. The second-order interpolation is used to transfer temperature at the gas-solid interface. Finally, the internal temperature distribution and the variation of thermal conductivity of gas-solid coupling are systematically studied. The main conclusions are as follows: (1) When the gas-phase scale is 60 nm, the boundary of the gas-phase region shows obvious discontinuous temperature distribution characteristics, and the temperature jump is obvious; when the transition of the temperature distribution in the solid-phase region is gradually uniform, the solid-phase scale changes from 2 nm to 5 nm.
(2) When the solid-phase scale is 5 nm, the temperature jump at the gas-solid interface decreases with the increase in the gas-phase scale. The temperature jump decreases sharply when the gas-phase scale is 10-70 nm and begins to change gently when the gasphase scale reaches 70 nm.
(3) The variation of thermal conductivity of gas-solid coupling with solid-phase scale is basically the same. The variation range of thermal conductivity is relatively small when the solid-phase scale is less than 1 nm. When the solid-phase scale is larger than 1 nm, the variation range of thermal conductivity begins to increase rapidly.
(4) The thermal conductivity of gas-solid coupling increases gradually with the increase in the gas-phase scale. The thermal conductivity increases greatly between 10 and 70 nm. When the thermal conductivity reaches about 70 nm, the increase in thermal conductivity begins to slow down gradually.  Data Availability Statement: The study did not report any data.

Conflicts of Interest:
The authors declare no conflict of interest.  Figure 11. Dependence of Coupled Thermal Conductivity on Solid-Gas Scale.

Conclusions
In this paper, the gas-solid coupling heat conduction model of nanoporous aerogels is established. Based on the theory of regional splitting, the gas-solid two-phase simulation is carried out by using the respective LBM. The second-order interpolation is used to transfer temperature at the gas-solid interface. Finally, the internal temperature distribution and the variation of thermal conductivity of gas-solid coupling are systematically studied. The main conclusions are as follows: (1) When the gas-phase scale is 60 nm, the boundary of the gas-phase region shows obvious discontinuous temperature distribution characteristics, and the temperature jump is obvious; when the transition of the temperature distribution in the solid-phase region is gradually uniform, the solid-phase scale changes from 2 nm to 5 nm.
(2) When the solid-phase scale is 5 nm, the temperature jump at the gas-solid interface decreases with the increase in the gas-phase scale. The temperature jump decreases sharply when the gas-phase scale is 10-70 nm and begins to change gently when the gas-phase scale reaches 70 nm.
(3) The variation of thermal conductivity of gas-solid coupling with solid-phase scale is basically the same. The variation range of thermal conductivity is relatively small when the solid-phase scale is less than 1 nm. When the solid-phase scale is larger than 1 nm, the variation range of thermal conductivity begins to increase rapidly.
(4) The thermal conductivity of gas-solid coupling increases gradually with the increase in the gas-phase scale. The thermal conductivity increases greatly between 10 and 70 nm. When the thermal conductivity reaches about 70 nm, the increase in thermal conductivity begins to slow down gradually.