Heat transfer and Entropy analysis of airflow around an airfoil subjected to an external magnetic field

In the present study, the thermal Lattice Boltzmann Technique is combined with the Ghost Fluid method to simulate the flow and heat transfer rate around a NACA 0015 airfoil in the presence of an external magnetic field source. It is tried to investigate the mutual effects of the Hartmann and Reynolds numbers as well as the magnetic field angle and the attack angle of the airfoil on the flow and heat transfer characteristics. Besides, the total entropy generation rate of the system was studied through its main components, i.e. entropy generation rate due to friction, heat transfer, and Magneto Hydrodynamics. Therefore, the tests were carried out for various Re and Ha numbers, plus different magnetic field angles and airfoil attack angles, and their influence on the active parameters which are Cd, Cl, and Nu, along with the entropy generation rate of the system,were recorded. The results revealed that with an increment of the Re number the Cd and Cl graphs drop, but the Nu value raises. Also, the total entropy generation rate of the system is at its maximum around 60   . fictitious equilibrium schemes, Interpolation schemes, Non-equilibrium extrapolation scheme, and the Ghost fluid schemes.


Introduction
A lot has changed since Qian et al. [1] first introduced the Lattice Boltzmann Method and it promptly became one of the most popular simulation tools. Its easy implementation, accurate results, and numerical efficiency were proved to be outstanding compared to other classical CFD techniques [2]. However, just like any other computational method,LBM suffers from various drawbacks and challenges [3][4].The first challenge that has been widely attended to, was to implement a suitable and efficient scheme to model the energy equation along with the Navier-Stokes equations. According to Guo and Shu [5], three different methods for the energy equation have gainedthe highest popularity so far. These schemes consist of the Multi-speed scheme, the Hybrid techniques, and the Double Distribution Method.
Multi-speed techniques apply different discrete velocities on the same lattice of the fluid flow, and so the particles can move to farther neighboring nodes. These methods are divided into low-order and highorder models [6][7]. However, according toChen et al. [8], some of the Low-order techniques suffer unwanted artifacts. The Hybrid methods treat the energy equation explicitly using Finite Volume or Finite difference schemes. Between the two, the latter gained more popularity [9][10][11] since it uses the same uniform grids as the fluid flow lattice. Although most of these Hybrid techniques have stability issues, recently Ahrar and Djavareshkian [12]proposed a TVD 1 FD-LBM 2 technique, which was successfully applied to High Ra problems. Also, the Double Distribution Function techniques are widely applied in the open literature [13][14][15]. In which the flow and energy equations utilize two different distribution functions. These schemes are considered as highly stable and accurate, numerically efficient, and easy-to-apply methods.
The other challenge that many researchers are dealing with, is the curved boundary treatments in the Lattice Boltzmann Method [16]. LBM uses a uniform grid in both flow and energy fields by definition and hence since the grid is not fitted to the obstacle surface, the curved boundaries are not fully captured by the technique.In order to resolve the predicament, several approaches are proposed in the last decade. Up to this day, 5 different main treatments are available in the open literature [17], namely: the Bounce-back 2 schemes, fictitious equilibrium schemes, Interpolation schemes, Non-equilibrium extrapolation scheme, and the Ghost fluid schemes.
The Bounce-back schemes which arevery popular due to their easy implementation in rectangular geometries are divided into two subdivisions. The nodal bounce-back scheme is applied to the node nearest to the physical boundary and the unknown distribution function is approximated by a distance fraction value [18]. On the other hand, the link bounce-back method specifies the mid-point of the link between the fluid node and the physical wall and similar to [19] the bouncing process occurs at this point. In either scheme, the curved wall is approximated by a staggering stair step and thus their accuracy is degraded to the firstorder. The Fictitious equilibrium scheme which was first proposed by Filippova and Hanel [20] was one of the first schemes that could capture the actual accurate curved wall geometry. The scheme applies the bounce-backed and fictitious distribution functions to estimate the post-collision distribution function and then the streaming process is carried out. This scheme also suffers severe instabilities in low-relaxation-time problems [21].
While the fictitious scheme also applies the same principles of the bounce back schemes, the Interpolation schemes which first came to the light by Bouzidi et al. [22]use a linear interpolation function to compute the unknown distribution function. The schemes are assumed to possess second-order accuracy in general but the accuracy can be enhanced by high-order interpolating functions. Unlike the previous schemes,the non-equilibrium extrapolation scheme carries out the collision step at the solid node. Guo et al. [23] proposed the scheme in 2002 and stated that the scheme benefits good computational stability as well as second-order accuracy.
In 2012 Tiwari and Vanka [24] introduced the Ghost Fluid Technique to Lattice Boltzmann method and claimed that their method can be applicable for all types of boundary conditions. In this scheme, at first, the ghost nodes, nodes that are located in the solid wall but have at least one link with the fluid nodes, are identified and then the equilibrium particle distribution functions at the ghost points are calculated from the extrapolated velocity and density values, and the non-equilibrium particle distribution functions are extrapolated from the fluid domain. The scheme is second-order accurate, applicable to both stationary and moving boundaries, and due to its comprehensiveness (It can also be utilized in classical CFD techniques like FVM 1 or FDM 2 [25][26]), it was able to gain wide popularity amongst the authors for different problems and various geometries [27].
In 2017 Mozafari-Shamsi et al. [28]applied the Ghost fluid Lattice Boltzmann Method (GF-LBM) to the natural convection of air between two concentric and eccentric cylinders. In addition to that, they simulated the sedimentation process of a cold cylindrical object in the fluid flow. In the same year, Chen et al. [29] applied a modified ghost fluid technique to reduce the spurious pressure oscillations on moving walls. According to their findings, this modified method was successful at decreasing the pressure oscillations while it was able to maintain the accuracy of common GF methods. A year after that, Ravindra et al. [30]examined the GF 3 technique on a cascaded lattice Boltzmann method to analyze the mixed convection heat transfer over a circular cylinder. They stated that the Drag and Lift forces as well as the heat transfer rate of the cylinder increase by an increment in Ri 4 number.
More recently Tong and Xia [31] explored the reliability of SRT 5 and MRT 6 Lattice Boltzmann Method with 3 different approaches for the curved boundary treatments. Their results revealed that the MRT technique could improve the numerical stability of the code as well as the accuracy of the results with the Ghost Fluid technique. Also,Balotaki et al. applied the ghost fluid technique along with double distribution function thermal LBM to simulate the natural convection of various nanofluids in a triangular cavity [32]. In this study, the surface Nusselt number and total entropy generation rate of the system were also observed. Saffarzadeh and Djavareshkian [33] simulated the flow over an airfoil using GF-LBM in 2020. To the best of our knowledge, the GF-LBM technique was usually tested against heat transfer problems and there is a lack of the validity and applicability of the technique for hydrodynamic cases in the open literature. Therefore, in this study, it is tried to present a comprehensive study of the GF-LBM technique for simulation of flow and the heat transfer rate of air over a NACA-0015 airfoil in the presence of an external magnetic field source with different inclination angles. Moreover, not only the mutual effect of the magnetic field angle and the airfoil angle of attack was observed but also the entropy generation rate of the system was examined for different cases of Re and Ha numbers.

Lattice Boltzmann Method
In the present study, the Navier Stokes equations are solved along with the energy equation to obtain the temperature field around the airfoil. The dimensionless equations are as follows: Re sin cos sin Re Re sin cos cos Re Here the third terms on the right-hand side of equations 2 and 3, are due to the Lorentz force induced by the external magnetic field source with an inclination angle of  and magnetic field intensity of 0 B . It is also worth mentioning that the induced magnetic field due to the magnetization of the fluid is considered to be very small with regards to the main magnetic field and the Buoyancy effects, as well as the Joule heating, are neglected in what follows. So, these dimensionless equations are to be solved to obtain the fluid flow and heat transfer results. In the present study, the well-known LBM technique is applied with the standard two-dimensional, nine-link (D2Q9) system for the temperature and flow fields and thus a brief review on the method is presented here but for a more detailed discussion, one can refer to [34]. The main discrete equations that are solved here, are as follows: In this regard, eq. (5) will recover the continuity and momentum equations (eq. 1-3) and eq. (6) models the energy equation (eq. 4). In these equations F representsany external force that can act on the flow domain (which in this case is the Lorentz Force), andthe equilibrium distribution functions eq i f and eq i g are calculated from: In these equations, c is equal to / xt  , with x  and t  being the space between two successive nodes and lattice time step, respectively. The weighing factors are also calculated from: Finally based on the Chapman-Enskog theory the mesoscopic relaxation times (   and   ) for the equations (5) and (6)

Ghost Fluid Method
In order to model the curved boundaries a bilinear interpolating ghost fluid scheme is applied. The main algorithm of the scheme is based on the study conducted by Tiwari andVanka [24]. In the first step, the ghost points must be identified. A ghost point is a node out of the fluid domain that has at least one link with one of the fluid nodes. Hence in the streaming process, there would be at least one unknown distribution function in the fluid nodes regarding this point. Secondly, an Image point must be found for each of the ghost points. An image point can also be described as a point on the wall-normal vector with the same distance of the ghost point with the wall and in the fluid domain. Then, the main flow parameters must be achieved at the image points by an interpolation of the values on neighbor nodes. In this regard a bilinear interpolation method is used: Here,  can be any of the flow main parameters such as T,  , u, and v. In the bilinear interpolation method, the four coefficients are the unknowns and must be calculated. To obtain the coefficients four different equations are required that can be written from the x, y, and  values at the four neighbor nodes.
In equation (14), the prime symbol represents the midpoints values, while i a can be calculated from: In the next step the calculated values for the image points will be extrapolated for the ghost nodes via a linear extrapolation technique: And finally, this distribution function will be applied in the streaming process to simulate the curve boundary no-slip condition.

Entropy Analysis Scheme
It is well known that any regular system can be subjected to several types of irreversibility depending on its situation and environment. In this study, three main types of entropy generation rates are considered, namely: The Frictional, Heat transfer, and MHD entropy generation rates.According to the local thermodynamic equilibrium of the linear transport theory [36]and regarding the dimensionless parameters for a 2D-Cartesian coordinate the entropy generation rates can be computed from: , which is assumed to be constant in this study and it is equal to 10 -4 . Also, the total entropy generation can be calculated from: Moreover, to observe the heat transfer and MHDentropy generations' importance with regards to the frictional irreversibility, the local Bejan number is introduced: This parameter can define the relative importance of heat transfer and MHD irreversibility in the total entropy generation rate of the system. Also, to study the overall dominant irreversibility mechanism of the system, the average Bejan number is computedfrom:

Geometry definitions and Validity Analysis
In the present study, a NACA 0015 airfoil is considered in a rectangular enclosure with a constant velocity and temperature profiles at the inlet. The outlet boundary condition is assumed to be a usual outflow condition that satisfies the continuity equation for the flow and the no-slip condition is valid for both horizontal walls. Moreover, the airfoil is stationary with a constant temperature boundary condition on the surface and the whole flow domain is subjected to a constant magnetic field source as can be seen in figure 2. In what follows the steady-state condition is considered and effects of the induced magnetic field, as well as the Joule heating, are neglected.
In order to test the reliability and accuracy of the results, various validation and mesh independence analyses have been carried out.The validation results are presented for flow past a 2d circular cylinder with Pr=0.71 and Re numbers of 10, 20, and 40. Table 1 indicates that the results of the Cd and Nu number of the present study are in excellent accordance with the published literature. Additionally, the mesh independence analysis was carried out for figure 2 geometry and 6 different mesh combinations with 367000, 480000, 607000, 750000, 1080000, and 1470000 grids were applied. As can be seen in figure 3 less than 0.1 percent of alteration was observed between the two last combinations and hence 1080000 grid-mesh was used in this study. Dennis and Chang [37] 1.522 _______ He and Doolen [38] 1.499 _______ Nieuwstadt and Keller [39] 1.550 _______ Sharma and Dhiman [40] _______ 3.237 Bharti et al. [41] _______ 3.282 Present study 1.573 3.124 Figure 2: the defined geometry of the airfoil and the external magnetic field. The airflow is simulated for the mentioned airfoil with the predefined geometry and conditions. At first, the fluid flow and heat transfer of the domain are solved in the absence of the external magnetic field. Figure 4 demonstrates the temperature contours of the flow around theairfoil at Re=200 and for different attack angles of 0, 6, 12, and 15. It can be seen that as the attack angle of the airfoil increases the wakes appear at the trailing edge of the airfoil. These wakes are released one after each other to produce a lowpressure domain behind the airfoil. Apparently, as the attack angle raises the pressure drag gains power over the frictional drag. On the other hand, the lift force is almost negligible at 0   due to the symmetry of the airfoil, but as the attack angle increases the lift force also augments. Besides, the same behavior is observed for the heat transfer rate and its representative Nusselt number. Figure 4: The temperature contours for the airfoil at Re =200 and Ha =0 with attack angles of 0, 6, 12, and 15. Figure 5presents the Cd, CL, and mean Nu graphs for various attack angles (  ) and three Re numbers of 120, 160, and 200. As can be seen, the Cd and CL graphs offer higher values for lower Re numbers, while the Nu graph shows a different manner. Since Nu number has a very close and direct relationship with the velocity magnitude, as the Re number rises, Nusselt number augments, correspondingly. It can also be noticed that the CL graph acts differently in low and high attack angles. Although in attack angles lower than 9 the lift coefficient of the highest Re number exceeds other cases, but as we approach the stalling point, it is clear that the altitude of the plots of the lower Re numbers, increases more drastically. Besides, it can be seen that almost in all cases there is a maximum point on the   Figure 6illustrates the Cd, Cl, and Nu graphs for various Ha numbers at Re=40, 80, and 120. In these plots, the attack angle of the airfoil and the magnetic field inclination angles are fixed at 6 and 45 degrees, respectively. It goes without saying that almost in all cases Cd and Nu increase with an increment of Re number. Moreover, it is seen that the drag coefficient rises with Ha augmentation. This can be due to the reduction of the velocity magnitude in the domain caused by the Laurentz forces and their braking natureon the flow. This phenomenon can be also noticed in Figure 6-c by a mild decrement of the Nu number. Meanwhile, Cl shows quite different manners in lower and higher Ha numbers. In Ha numbers less than 10 it is seen that the lift coefficient is higher in the lower Re numbers, however, in larger Ha numbers the effect of the stream function strength gives rise to the Cl magnitude. Moreover, although they are different parameters it seems that the magnetic field and attack angles can have a mutual effect on the flow behavior. To study their influence on the aerodynamics and heat transfer rate of the airfoil Figure 7 presents Cd, Cl, and Nu number of the airfoil in different magnetic field angles and attack angles. It is seen that almost in all magnetic field angles the Cd and Cl rise with an increment of attack angle. Also, with an increment of the magnetic field angle Cd and Nu graphs reach a maximum point and then their values reduced. This maximum can happen because of better air circulation and proves that the induced magnetic field loses dominancy at these angles.

Entropy generation analysis of the system
In this section, the entropy generation of the system is studied for different cases. Figure 8 presents the entropy generation of the flow in the absence of the external magnetic field. As the Re number increases, the heat transfer and frictional entropy generation augment due to the higher velocity and temperature gradients in the domain. On the other hand, due to the nature of forced convective heat transfer, it is obvious that the heat transfer irreversibility is much stronger than the frictional entropy generation.
Besides, there is a maximum point on all graphs at 12   , which is the same attack angle representing the maximum points of graphs in Figure 5. At this attack angle, the highest rate of stream function is achieved and hence the rate of heat transfer as well as its corresponding entropy generation rate occurs. Moreover, to investigate the magnetic field inclination angle influence on the entropy generationrate of the system Figure 9 presents different components of the entropy generation rate of the system as well as the total rate of irreversibility production for different attack angles and Re=100 and Ha=15. It can be seen from Figure 9-a and b that the rate of entropy generation rises with an increment of the attack angle. However, the magnetic field angle also plays an important role in the entropy generation of the system. It is seen that almost for any of the attack angles, the minimum entropy generation due to the heat transfer and friction occurs at 0   and the maximum point is somewhere around 45   . This can be related to the fact that at this point the air circulation is much better, hence the velocity of the flow, and consequently it gives rise to the entropy generation rate of the system.
It is also worth mentioning that unlike these two graphs, the MHD entropy generation rate Figure 9-c was not strongly dependenton the attack angle. It can be derived from equation (20) that the MHD irreversibility only is affected by the Ha number and variation of two flow velocity components.Thus, even though the latter is influenced by the attack angle, but it is obvious that the effect of the former is more pronounced in the final result. Finally, the total entropy generation of the system Figure 9-d is strongly under the influence of its heat transfer component and it is mildly affected by the frictional and MHD entropy generations. Figure 9: the heat transfer entropy generation (a), frictional entropy generation (b), MHD entropy generation (c), and total entropy generation rate of the system (d) for different magnetic field angles at Re=100 and Ha=15.

Conclusion
In this study the influence of an external magnetic field source on the heat transfer and entropy generation rate of a NACA 0015 airfoil has been studied, numerically. The well-known LBM technique was incorporated with the Ghost Fluid method to simulate the flow around a curved boundary. Also, the magnetic field intensity and direction as well as the attack angle of the airfoil were altered and their influence on the Cd, Cl, and Nu was recorded. The main findings of the present research are as follows:  As the Re number increases, the Cd and Cl coefficients are reduced, but the Nusselt number and the heat transfer rate of the airfoil increase in the absence of an external magnetic field.  With an increment of the Ha number and due to the Lorentz force the flow stream function drops and this will lead to a decrement of the mean Nu number.  The study of the mutual effects of the attack angle and magnetic field angle proves that generally, the maximum values of Cd and Nu occur at 0   , however, the Cl graph has its maximum value at 45   .
 The heat transfer and frictional entropy generation rates of the system increase drastically with an augmentation of the attack angle.  Finally, it was observed that for different magnetic field angles the maximum entropy generation rate is found at around 60   , where the highest rate of stream function is achieved.