Application of Response Surface Methodology and Artificial Neural Network to Optimize the Curved Trapezoidal Winglet Geometry for Enhancing the Performance of a Fin-and-Tube Heat Exchanger

Abstract

The present work aims at optimizing the geometry of curved trapezoidal winglets to enhance heat transfer rates (expressed as Colburn factor, j) and minimize pressure losses (expressed as friction factor, f). A fin-and-tube heat exchanger was analyzed with winglets mounted on the alternate tube and on either side of the fins. Multi-objective optimization was performed using the genetic algorithm (GA) to maximize j and minimize f. Two surrogate models, viz. response surface methodology (RSM) and artificial neural network (ANN), were considered as inputs to GA. To reduce the number of runs, a sensitivity analysis was first performed to select the most influential geometrical parameters for optimization. The values of j and f in the design of the experiments table were computed using CFD. The Pareto front points elucidated a significant improvement compared with the reference model along with a broad choice for the designers, not only for the design condition but also for the off-design inlet condition.

1. Introduction

Fin-and-tube heat exchangers (FTHEs) find numerous applications such as in automobile radiators, air conditioning, petrochemical industries, and many more. Conventional heat exchangers suffer a major drawback of high space requirements and low airside heat transfer coefficient. Several heat transfer enhancement techniques, viz. the use of nanoparticles, corrugated surfaces, elliptical tubes, wavy fin, and vortex generators in the flow field, have been reported in the last few years to overcome these setbacks and increase the overall thermo-hydraulic performance of heat exchangers.

Recently, the application of vortex generators (VGs) to increase the heat transfer coefficient on the airside has become popular to enhance HEs’ efficiency [1,2]—these VGs produce intense longitudinal vortices without many frictional losses [3]. Actually, VGs installed on the airside increase heat transfer by redirecting the main fluid flow toward the recirculation zone. However, the increase in heat transfer is accompanied by enhanced pressure losses and this sets the foundation for the optimization of vortex generator geometry (in our case, curved trapezoidal winglets). VGs can be arranged in common flow-up and flow-down configurations; the former is known to be more efficient in accelerating the flow and delaying flow separation [3,4,5,6,7,8]. The tubes can be arranged in an inline or staggered array arrangement in an FTHE [9,10]. The shape of the tubes in an FTHE is generally circular, but elliptical tubes have also been reported in the open literature [5].

The geometry of the VG plays a major role in heat transfer by creating additional turbulence in the flowing fluid downstream. A comparative study was performed by [11] between curved rectangular winglets, sinusoidal wavy rectangular winglets, and rectangular winglets in an FTHE, and the authors reported that wavy rectangular VGs provided the highest heat transfer enhancement. In an FTHE, the fin pitch, longitudinal tube pitch, transverse tube pitch, length of vortex generator, the height of vortex generator, and attack angle of the delta winglet type VG were found to be the most important geometric entities in the view of thermal performance [12].

Apart from the straight profiles of the winglet, curved profiles have also gained popularity [13,14]. An experimental investigation was reported on the comparison of the shapes of the VG, viz. rectangular, delta, and trapezoidal with straight as well as curved profiles, by Zhou et al. [15] and via numerical simulation by Lu et al. [16] and Esmailzadeh et al. [17]. They reported that curved profiles are more promising to increase the thermal performance of the HE. In addition, a similar study was reported with punched holes near the base in VGs, wherein the pressure losses were reduced with a punched hole configuration [18,19,20]. The effect of the position and geometry of a curved VG in an FTHE was studied by Oh et al. [21]. A CFD study on an FTHE having a curved rectangular VG and staggered tube arrangement was reported by Gong et al. [19], with an aim of predicting the optimum design parameters, viz. circumferential position, radial position, base length, edge height, and fin spacing. Apart from a surface-mounted VG, some modified fin structures such as brazed plate HE [22], double-pipe HE [23], finned tube [24], and microchannel heat sink [25] have also been reported. One such surface modification in the form of arc-shaped baffles, studied and optimized by Promvenge et al. [26], reported that a 90° angle of arc-shaped baffles and a pitch ratio of 8 offered the highest heat transfer enhancement. The impact of tube diameter and fin pitch on the Nusselt number (Nu) in an FTHE with a punched delta winglet VG and punched rectangular VG was investigated by Wu et al. [27] and Tian et al. [28]. The experimental results revealed that the fin pitch had a more prominent effect on Nu compared to the tube diameter.

Apart from analyzing the impact of variations in an individual design parameter, multi-objective optimization studies have also been reported in the open literature. The latter is advantageous to achieve the global optimized value(s) of the design parameters to enhance the overall thermal performance of the HE. Different optimization techniques have been proposed and used by various researchers to achieve the goal, viz. RSM (response surface methodology) [29], ANN (artificial neural network), and GA (genetic algorithm). An optimization study utilizing RSM was reported by [30] for the optimization of the Ranque-Hilsch vortex tube refrigerator and found that the number of intake nozzles and orifice diameter largely influenced the rate of heat transfer. Hosseinirad et al. [31] performed a CFD investigation for the analysis of wavy VG and used RSM for its optimization. Kanaris et al. [32] employed RSM and developed a correlation between the design parameters and response variables based on the CFD results, and later on predicted the optimal design for the plate HE. Khaljzadeh et al. [33] performed a similar study using RSM on vertical ground HE.

A comparative study between the RSM- and ANN-based optimization techniques applied in a V-down rectangular channel was reported by Chamoli et al. [34]. It was further stated that RSM is more accurate in predicting, with a maximum deviation of 5%. Apart from RSM, several researchers performed combined analyses using CFD, RSM, and GA. For instance, Lemouedda et al. [35] optimized delta winglet VGs in a fin-and-tube HE for an optimum VG shape and attack angle. A similar investigation was performed by Sun et al. [36] but for elliptical tubes instead of circular ones. Gholap et al. [37] optimized the energy consumption and material cost in an HE for refrigerators. The effect of fin spacing, the ratio of longitudinal to transverse tube pitch, and the attack angle of the delta winglet in an FTHE were studied by Wu et al. [38] using RSM, and optimized using GA. Similar optimization techniques have been used by several authors, for instance, for the analysis and optimization of a novel helical coil tube in a tube HE by Liu et al. [39], and in a conical strip VG by Zheng et al. [40]. A comparative study between RSM and direct optimization (DO) was performed by Salviano et al. [41] to optimize the VG angle and position in an FTHE and they also studied the influence of VG roll angle on heat transfer. They reported that the DO approach delivered better results than the RSM approach.

The artificial neural network (ANN) approach has gained immense popularity in the last few decades due to its capability to tackle non-linearity with great accuracy [42,43]. The work of Haykin [44] and Eeckman [45] provides a detailed history and mathematical background of neural networks. Berber et al. [46] used ANN for the prediction of heat transfer in a circular tube. A double-pipe HE with delta-wing-type inserts was studied by Kocyigit et al. [47] and Khan et al. [48] and they used ANN to validate the results with experimental data. Several authors used ANN as a surrogate model and GA for optimization. For instance, CFD, ANN, and GA have been clubbed together for the optimization of wavy fin-and-tube HE by Li et al. [6], who reported that an increase in corrugation height and angle of attack have a positive impact on the heat transfer. Additionally, a wavy fin height of 1.8 mm and an attack angle of 60° produced a maximum heat transfer effect. In a similar study, Abdollahi et al. [3] optimized the winglet geometry. Khan et al. [49] studied the trapezoidal and rectangular VG, attack, and attach angle in a plate-fin channel and concluded that an attack angle between 45 and 60° and an attached angle of 60° are the most suitable for maximum heat transfer and minimum pressure losses. Song et al. [8] studied a wavy fin with a circular tube HE equipped with a curved VG, and a maximum enhancement in heat transfer of 30.6% was achieved. Instead of surface modifications, some researchers opted for an active heat transfer enhancement method, introduced nanofluid for increasing the heat transfer performance [50], and achieved a maximum enhancement of 82%. Ke et al. [51] considered winglet attack angle for optimization to enhance the performance of the HE.

From the literature, it can be confirmed that numerous experimental and numerical studies have been performed to analyze and ultimately enhance the heat transfer performance of heat exchangers equipped with vortex generators. The experimental and numerical studies on curved winglets [52,53] are limited and the scope for multi-objective optimization is vast. Thus, the present focus is on a multi-objective genetic algorithm (MOGA) with an intention of increasing the heat transfer capability and minimizing the pressure losses, and thus generating a set of optimal data points (Pareto fronts) that best satisfies the objective.

Although RSM/ANN has been used in some of the multi-objective studies, the present work incorporates both of these approaches (separately) for establishing a relation between design parameters and objective function and then applying MOGA to compute Pareto fronts. Although Sarangi et al. [54] performed multi-objective optimization of the curved trapezoidal winglets, they considered only three parameters for optimization, viz. trailing edge height, arc radius, and the angle subtended on the tube center. Apart from these parameters, the present work accounts for the winglet front edge height and the eccentricity of winglets for analysis. Furthermore, the most sensitive parameters are considered for optimization.

To summarize, the theme of the present study is to enhance the overall thermal performance of the FTHE by optimizing the curved trapezoidal winglet geometry. To achieve this, the work is planned as follows: first, the flow domain is modeled using curved trapezoidal winglets on alternate fin walls placed on either side. Second, different levels of mesh are generated and tested to ensure the results are grid-independent, followed by a validation study to justify the solver settings and boundary conditions used. Third, before proceeding with optimization, the critical parameters having a significant effect on the performance of the HE are identified. Fourth, based on these parameters, the design of the experiment (DoE) is prepared using the factorial design method. Fifth, we make use of response surface methodology as well as the Bayesian-regularized ANN surrogate model as the two separate surrogate models. Sixth, a multi-objective optimization study is performed using GA to maximize the heat transfer and minimize the pressure losses, and two separate Pareto fronts are generated using RSM and ANN surrogate models. Finally, a few Pareto fronts, selected from both of the surrogate models, are selected and compared with the CFD results to measure the relative error and explore the flow physics of the optimized models.

2. Model Description

2.1. Geometry of Fin-and-Tube Heat Exchanger

In the present work, we used the HE geometry proposed by Joardar and Jacobi [55] as the reference model—the same was also used for validation. Different experiments were carried out by Joardar and Jacobi [55] on HEs, having seven inline rows of tubes. The first HE did not incorporate any winglets and this model was taken as the baseline model (BLM) here. They also considered one and three pairs of the delta winglets mounted on one and three (alternate) tubes, respectively, and these models were compared to the optimized models. Regarding the geometrical details, the tubes have a diameter (D) of 10.67 mm, transverse pitch (P_t), and longitudinal pitch (P_l) of 25.4 mm each. The length (L) of the FTHE was 177.8 mm, the fin pitch (F_p) was 3.63 mm, and the fin thickness (t) and winglet thickness were 0.18 mm each. Instead of delta winglets—as taken by Joardar and Jacobi [55]—we considered curved trapezoidal winglets mounted alternately on either side of the fin surface as well as on the alternate tubes. A detailed view of the geometry is presented in Figure 1.

Figure 1a represents the geometry of the entire model. Owing to its symmetrical nature, we considered a small part of the entire domain (marked by the enclosure). Due to the repetitive nature of the geometry about the plane containing the fin surface (cf. Figure 1b), a periodic boundary condition was applied to these plates (with half the thickness of the original fin plate). The implementation of the two approaches largely reduced the size of the geometry, and hence the computational overheads. As mentioned earlier, the winglets were mounted on either side of the fin plate (cf. Figure 1b) and the alternate tubes (cf. Figure 1c).

2.2. Equations Used for Calculating the Performance of HE

The performance of an HE is expressed in a dimensionless form and computed based upon several parameters, in accordance with Joardar and Jacobi [55], as follows:

$$\mathrm{Total} \mathrm{heat} \mathrm{transfer}, Q=\dot{m}{C}_p\left({\overline{T}}_o-{\overline{T}}_{in}\right)$$

(1)

$$\mathrm{Logarithmic} \mathrm{mean} \mathrm{temperature} \mathrm{difference}, LMTD=\frac{\left({\overline{T}}_{in}-T_w\right)-\left({\overline{T}}_o-T_w\right)}{\ln \left(\frac{{\overline{T}}_{in}-T_w}{{\overline{T}}_o-T_w}\right)}$$

(2)

$$\mathrm{Minimum} \mathrm{free} \mathrm{flow} \mathrm{area}, A_{min}=\left(B-D\right)\times Fp$$

(3)

$$\mathrm{Total} \mathrm{heat} \mathrm{transfer} \mathrm{surface} \mathrm{area}, A_T=2BL+7\pi DFp$$

(4)

$$\mathrm{Heat} \mathrm{transfer} \mathrm{coefficient}, h=\frac{Q}{A_{T}LMTD}$$

(5)

$$\mathrm{Hydraulic} \mathrm{diameter}, D_h=\frac{4A_{min}L}{A_T}$$

(6)

$$\mathrm{Reynolds} \mathrm{number}, Re=\frac{\rho U_{in}{D}_h}{\mu }$$

(7)

$$\mathrm{Pressure} \mathrm{drop}, ∆p=p_{in}-p_{out}$$

(8)

$$\mathrm{Colburn} \mathrm{factor}, j=\frac{h}{\rho U_{in}{C}_p}P{r}^{2/3}$$

(9)

$$\mathrm{Friction} \mathrm{factor}, f=\frac{∆p}{0.5\rho U_{in}^2}\frac{A_{min}}{A_T}$$

(10)

$$\mathrm{Enhancement} \mathrm{factor}, \eta =\left(\frac{j}{j_o}\right)\times {\left(\frac{f}{f_o}\right)}^{\left(\raisebox{1ex}{$-1$}\!\left/ \!\raisebox{-1ex}{$3$}\right.\right)}$$

(11)

Here, $\dot{m}$ is the mass flow rate in kg/s, $C_p$ is the specific heat capacity in J/kgK, ${\overline{T}}_o$ and ${\overline{T}}_{in}$ are the area-weighted average temperature at the outlet and inlet, respectively, in Kelvin (K), $T_w$ is the wall temperature, $\rho$ is the density of air in kg/m³, $U_{in}$ is the inlet velocity in m/s, $\mu$ is the kinematic viscosity in Pa·s, $p_{in}$ and $p_{out}$ are the pressure at the inlet and outlet, respectively, in Pascal (Pa), and $Pr$ is the Prandtl number defined as the ratio of kinematic viscosity (μ) and thermal diffusivity (α). B is the width of the computational domain and is equal to 25.4 mm. The Colburn factor (j) was used to calculate the heat transfer and the friction factor (f) was used to represent the pressure drop. The enhancement factor, given by Equation (11), was adopted from Webb and Kim [56], where j and f are calculated based on the selected Reynolds number (Re), and j_o and f_o are reference values taken from the baseline model under similar operating conditions. Using these reference values in the equation makes it a non-dimensional parameter, useful for calculating the rating.

2.3. Description of Numerical Simulation

Across the last few decades, computational fluid dynamics (CFDs) have gained tremendous popularity and have been widely used as a robust tool for thermal and fluid flow analysis [57]. CFD is a well-proven numerical technique and the best alternative to experimentations due to significantly fewer capital investments, higher accuracy, and lesser time requirement. If performed carefully, the CFD calculations are as good as the experimental results. Furthermore, with CFD, we have the advantage to visualize and extract the flow field details at any time. It also provides the flexibility to make changes in the geometry during design optimization with no additional cost factor involved, other than the power consumption by the hardware.

The present study makes use of Menter’s transition shear stress transport (SST) model—a four-equation model—together with curvature correction for modeling turbulence. This transition model utilizes the two equations of the SST turbulence model along with two additional transport equations, viz. gamma (γ) for determining intermittency and R_eθ for determining laminar–turbulent transition process determination [58].

The following governing equations are solved for the fluid zone:

Continuity equation:

$$\frac{\partial }{\partial x_i}\left(\rho u_i\right)=0; i = 1,2 \mathrm{and} 3$$

(12)

Momentum equation:

$$\rho \left(u_i\frac{\partial u_j}{\partial x_i}\right)=-\frac{\partial P}{\partial x_j}+\frac{\partial }{\partial x_i}\left[\left(\mu +{\mu }_t\right)\frac{\partial u_j}{\partial x_i}\right]; i, j = 1,2 \mathrm{and} 3$$

(13)

Energy equation:

$$\frac{\partial }{\partial x_i}\left[u_i\left(\rho E+P\right)\right] = \frac{\partial }{\partial x_i}\left[{\lambda }_{eff}\frac{\partial T}{\partial x_i}\right]$$

(14)

Here, E signifies the total energy, and the effective thermal conductivity is given by ${\lambda }_{eff}$ = $\lambda$ + ${\lambda }_t$, where ${\lambda }_t$ is the turbulent thermal conductivity.

The transition SST model uses the standard SST model as the base, given in terms of the turbulent kinetic energy (k) and specific rate of turbulence dissipation (ω), as follows:

$$\frac{\partial \left(\rho u_{i}k\right)}{\partial x_i}={\widehat{P}}_k-{\widehat{D}}_k+\frac{\partial }{\partial x_i}\left[\left(\mu +{\sigma }_k{\mu }_t\right)\frac{\partial k}{\partial x_i}\right]$$

(15)

$$\frac{\partial \left(\mathit{\rho }{u}_i\omega \right)}{\partial x_i}=P_{\omega }-D_{\omega }+\frac{\partial }{\partial x_i}\left[\left(\mu +{\sigma }_{\omega }{\mu }_t\right)\frac{\partial \omega }{\partial x_i}\right]+2\left(1-F_1\right).\frac{\rho {\sigma }_{\omega 2}}{\omega }.\frac{\partial k}{\partial x_i}.\frac{\partial \omega }{\partial x_i}$$

(16)

The two transport equations relating to the intermittency and transition momentum thickness Reynolds number, coupled with the original SST model, i.e., Equations (15) and (16), respectively, are given as follows:

$$\frac{\partial \left(\mathit{\rho }{u}_j\gamma \right)}{\partial x_j}=P_{\gamma }-E_{\gamma }+\frac{\partial }{\partial x_j}\left[\left(\mu +\frac{{\mu }_t}{{\sigma }_f}\right)\frac{\partial \gamma }{\partial x_j}\right]$$

(17)

$$\frac{\partial \left(\mathit{\rho }{u}_j{\widehat{R}}_{e\theta t}\right)}{\partial x_j}=P_{\theta t}+\frac{\partial }{\partial x_j}\left[{\sigma }_{\theta t}\left(\mu +{\mu }_t\right)\frac{\partial {\widehat{R}}_{e\theta t}}{\partial x_j}\right]$$

(18)

The energy equation in the solid regime reads as follows:

$${\gamma }_{ij}\frac{\partial }{\partial x_i}\frac{\partial T}{\partial x_i}=0$$

(19)

Here, ${\gamma }_{ij}$refers to the conductivity matrix.

For further details on the production term (P_k), destruction term ($D_k),$blending function (F₁), the turbulent viscosity (µ_t), and all of the other terms, including the details of the model constants, we referred to [59,60].

2.4. Discretization of the Solution Domain and Solver Settings

The first step in CFD is to discretize the entire domain into small control volumes (cells) over which governing equations are solved to generate a flow field. In the present work, we used hexahedral cells to discretize both fluids as well as solid domains. Ansys workbench was used to generate the Cartesian mesh over the entire volume and prism layers were provided over the cold tube walls to accurately capture the temperature gradients. The Cartesian cells were orthogonal, and hence, of high quality—they quickly filled up the entire domain and were trimmed off to conform to the solid boundaries. Furthermore, these cells are less diffusive in nature compared to the tetrahedral cells. Local grid refinement was conducted on the winglet walls using the octree method. Figure 2 shows the mesh details over the fins, cooling tubes, and VGs. The fin walls were made separately in ICEM CFD and appended as a case file in FLUENT, and interfaces were made between the solid and fluid domains manually so that conjugate heat transfer could take place.

Next, we performed a detailed study on the dependency of grid count on the performance parameters of the HE. Seven different levels of mesh were considered for analysis, and the percentage deviation from the experimental results is shown in

Table 1. Effect of cell counts on the performance, and deviation from the experimental data.

Element Size (mm)	Mesh Count (in Millions)	Re (-)	h (CFD) (W/m2K)	h (Exp. [55]) (W/m2K)	Variation w.r.t. BLM (%)
0.55	0.133	659.74	34.716	36.593	5.129
0.50	0.176		34.817		4.853
0.45	0.246		35.444		3.139
0.40	0.332		35.551		2.847
0.35	0.470		35.552		2.844
0.30	0.729		35.556		2.833
0.25	1.166		35.559		2.828

. It becomes apparent that after the fourth level of mesh, the deviation is less than 3%—the closest agreement is observed with the finer mesh consisting of more than 1 million cells. To obtain a good balance between the accuracy and the computational time, the grid topology corresponding to 332 K cells was selected for all future work.

Air was taken as a working medium with a density of 1.225 kg/m³ and viscosity of 1.885 × 10⁻⁵ kg/ms. At the inlet, the velocity inlet boundary condition was used, and a uniform velocity profile was prescribed with a mean velocity of 2.14 m/s, corresponding to a Reynolds number of Re = 634.5 at a temperature of 310.6 K. On the outlet plane, a pressure outlet boundary condition was used with the gauge pressure amounting to zero, and for backflow, a normal to boundary method was quoted. The temperature of the cold tubes was kept constant at 291.77 K. For aiding the conjugate heat transfer between fluid, fin walls, and winglets, the coupled boundary condition was implied. A no-slip boundary condition (u⸱$\widehat{n}$ = 0) was applied for the tube walls, fin surface, and winglets.

Referring to the schematic diagram of the HE (Figure 1), and as discussed earlier, a symmetry boundary condition was applied to the fluid zone over the plane described in Figure 1a—no flow takes place across these surfaces. Owing to the repeatability across the fin walls (cf. Figure 1b), a periodic boundary condition was applied. The application of these boundary conditions together reduces the size of the computational domain significantly with no loss in accuracy. The SIMPLEC (semi-implicit method for pressure-linked equation consistent) method was applied for the pressure velocity coupling. Standard discretization was implied for pressure. The QUICK scheme was implied for the discretization of turbulent kinetic energy, intermittency, momentum, specific dissipation rate, and momentum thickness Re and energy equations [61,62].

2.5. Validation and Methodology

The numerical model was validated against the experimental data of Joardar and Jacobi [55], consisting of a 3 VG pair array (cf. Figure 1, p. 1158). Figure 3 shows the variations in the heat transfer coefficient with different Re against the experimental data. A close agreement exists between the numerical calculations and the reference data—a maximum deviation of 8.94% (only at the highest Re) and a minimum deviation of 0.46% are observed.

2.6. Selection Criteria of Geometric Entities for Optimization

The effect of radius, the angle subtended by the winglet on tubes, the winglet leading as well as trailing edge heights, and eccentricity provided to the winglets upon heat transfer and pressure drop were studied. Here, one parameter is varied at a time while imposing three different values (to all of the other geometrical entities), namely, lower, middle, and upper bounds (also, cf. Table 2), as shown in Figure 4. The conclusive results indicate the following:

• With an increase in radius (r), both j and f values increase for the upper bound values of other parameters and this change is significant. For the middle bound values, the change in both j and f values is mild. At the lower bound, the j as well as f values decrease moderately until r = 9.5 mm and then increase with a further increase in r. This indicates the prominent impact of r at the upper bound values.
• The j and f values both increase with an increase in θ at the upper bound until θ = 60°, followed by nearly no change in the j value while the f value still increases. A similar trend is seen for the j value at the middle bound—initially, it increases until θ = 60° and then decreases by a small amount, while the f value increases consistently. At the lower bound, the j value decreases until θ = 45° and then increases mildly with θ, whereas the f value first undergoes reduction up until θ = 45° and thereafter it is nearly a constant.
• With an increase in the leading edge height (h₁), at the upper bound both the j and f values increase. At the middle bound, the j values first increase until h₁ = 2 mm and then remain constant while the f value continues to increase mildly. At the lower bound, the j value undergoes a more significant reduction than the f values do.
• The values of j and f increase with an increase in trailing edge height (h₂) at the upper bound and middle bound while this increase is significant at the upper bound. At the lower bound, the variations are quite small.

The effect of eccentricity on the j and f values is as follows: At the upper bound, the variations in the j values resemble a sinusoidal pattern of very low amplitude, while the f values elucidate a mild but continuous increment. At the middle bound for the j value, there is a slight increase up to e = 5 mm, followed by a decreasing pattern, and again a mild increase beyond e = 7 mm. The f values increase slightly with the eccentricity at the upper bound, while beyond e = 7 mm, a slight reduction is observed. However, in contrast to the upper bound results, at the lower bound, the f value (although small) reduces continuously with an increase in eccentricity values.

Figure 4. Variation in j and f values with changes made in the design parameters.

Figure 4. Variation in j and f values with changes made in the design parameters.

Table 2. Settings for the DoE.

Variable	Quantity	Lower Bound	Upper Bound
r	Winglet radius (in mm)	8	11
$\theta$	The angle subtended by winglet (in degrees)	30	75
h1	Winglet leading edge height (in mm)	0.7	3.3
h2	Winglet trailing edge height (in mm)	0.7	3.3

Table 2. Settings for the DoE.

Variable	Quantity	Lower Bound	Upper Bound
r	Winglet radius (in mm)	8	11
$\theta$	The angle subtended by winglet (in degrees)	30	75
h1	Winglet leading edge height (in mm)	0.7	3.3
h2	Winglet trailing edge height (in mm)	0.7	3.3

A general trend shows that with an increase in the values of the design parameters, the heat transfer and pressure drop both increase, particularly at the upper bounds. This is due to the fact that for the upper bound values of the geometrical parameters, a large amount of fluid from the mainstream is directed toward the wake region and helps in improving the mixing of the fluid. At the lower bound values, there is a slight decrease in heat transfer since the smaller geometrical values may not be able to influence the flow field to the extent of the upper bound. Interestingly, with the middle bound option, the performance undergoes a change after a certain value (as discussed earlier in detail).

The conclusive results indicate that the eccentricity does not contribute significantly to the heat transfer process, whereas the other geometrical entities demonstrate a remarkable effect on the performance of the HE (particularly at the upper bounds and to a certain value in the middle bounds). Hence, we consider four geometrical entities, viz. the winglet radius, the angle subtended by the winglet, and the height of the leading edge and trailing edge of the winglets. The range of values considered for optimization has been represented in Table 2.

2.7. Design of Experiments

Numerous researchers have performed optimization studies taking into consideration the effect of single design parameters while keeping the others fixed. This is a time-consuming and local approach, and the results obtained may not reach a globally optimized value. Therefore, the most widely accepted and adopted technique is to design the experiments in accordance with some standard methods. This approach ensures the participation of all of the chosen design parameters, and hence, is interactive in nature. The values of the dependent variables for each combination of independent variables are calculated by making use of either experimentations, numerical simulations, or (semi-)empirical relations. The next step is to fit the surface over these data to generate the equations, which are then subjected to optimization; hence, this step (to generate equation(s) prior to optimization) is often called a surrogate method. The accuracy of the results is highly dependent on the way that the experiment is designed and calculated, the method used for surface fitting, and the approach used for optimization.

Since the design space is large, a three-level factorial design—that divides each factor into three levels, viz. lower, middle, and upper levels—was used in the present study to analyze the non-linear interaction of the design variables. Considering the design parameters, a sample space consisting of 3⁴ = 81 design points was utilized to capture the quadratic relationship for each design variable. As mentioned earlier, four design parameters were chosen for optimization, and the extent of their variation is listed in Table 2. The chosen scheme fills the entire design space uniformly. The values of the Colburn factor (j) and friction factor (f) computed using CFD are listed for all of the design points in

Table 3. Input data for surrogate model.

No.	R (mm)	θ (deg.)	h1	h2	j (×103)	f (×103)
1	9.5	52.5	0.7	3.3	18.553	84.809
2	9.5	52.5	3.3	3.3	18.740	97.484
3	9.5	30	3.3	0.7	15.238	64.323
4	8	75	0.7	3.3	17.495	71.745
5	11	30	2	3.3	17.312	78.545
6	11	75	3.3	3.3	23.725	228.332
7	8	52.5	3.3	3.3	16.322	77.897
8	9.5	30	0.7	2	15.698	61.919
9	8	52.5	2	2	16.298	68.235
…	…	…	…	…	…	…
73	9.5	52.5	2	0.7	15.835	64.425
74	11	52.5	0.7	3.3	20.632	102.554
75	9.5	75	2	2	18.382	82.700
76	11	75	3.3	2	21.995	138.715
77	9.5	75	3.3	2	20.035	100.050
78	8	30	2	2	15.443	62.617
79	8	75	2	3.3	18.194	79.229
80	11	30	3.3	0.7	15.542	67.403
81	11	75	2	0.7	17.012	74.989

—since there are a large number of runs, only the first nine and last nine points are shown.

3. Results and Discussion

3.1. Surrogate Modeling with RSM

The RSM is an accurate method for the prediction of the relationship between input, output, and optimization of the design. It involves a set of mathematical and statistical tools that are used to model and analyze the relationship between input variables and output responses. The methodology involves a variety of experiments that are designed to evaluate the effects of different input variables on the output responses of interest and to identify the optimal combination of input variables that produce the desired response. RSM typically involves the use of statistical models such as regression models and analysis of variance (ANOVA) to analyze the data and identify the significant factors that affect the response variables.

The responses are graphically represented in the three-dimensional space via a surface plot. The response variable R for the second-order model is expressed as the sum of a constant term, a linear term (of input variables), the interactive terms (between the two variables), quadric terms, and an error term. It is represented as follows:

$$R={\alpha }_o+{\displaystyle \sum }_{n=1}^k{\alpha }_i{Y}_i+{\displaystyle \sum }_{n=1}^k{\alpha }_{ii}{Y}_i^2+{\displaystyle \sum }_{i<j}{\alpha }_{ij}{Y}_i{Y}_j+\epsilon$$

(20)

where R is the response variable, and ${\alpha }_0$, ${\alpha }_i$,${\alpha }_{ii}$, and ${\alpha }_{ij}$ are the regression coefficients for constant, linear, quadratic, and interactive terms, respectively. $Y_i$ and $Y_j$ are the independent variables and $\epsilon$ is the error term [29].

The main effect plot shows the impact of individual design variables on the objective function and helps to identify the most influential parameters. The effect of design parameters r, θ, h₁, and h₂ on j and f, respectively, has been shown in Figure 5. The slope and length of the curves signify the impact of the entity in concern. The direction represents the positive or negative influence on the performance [29]. With an increase in r, θ, and h₂, the heat transfer increases sharply, while with an increase in h₁, the heat transfer first increases and then flattens. The f values increase with an increase in design parameters and the trend is, somewhat, consistent for all of the independent variables.

Figure 6 and Figure 7 show the standardized effect of variables and their interactions on heat transfer and pressure drop, respectively, via the Pareto chart. Here, the vertical line (in red) corresponds to a 95% confidence level and the horizontal (bar) lines that cross this vertical line signify the most important terms (be it an individual term or interactive terms)—the length of the bar extending to the right from the vertical line represents its impact. From Figure 6 and the p-values in

Table 4. The p-values for j.

Term	p-Value
r	0.0000
θ	0.0000
h1	0.0000
h2	0.0000
r2	0.1296
rθ	0.0000
rh1	0.0002
rh2	0.0000
θ2	0.0038
θh1	0.0000
θh2	0.0000
h12	0.0150
h1h2	0.0000
h22	0.3526

R-squared = 97.256%. Bold numerals indicate significant factor at 95% confidence level.

, twelve significant terms at a 95% confidence level are observed. Interestingly, θ and h₂ are the most significant parameters, followed by r, and there is a comparatively lesser impact of h₁ on heat transfer—all four of them represent the individual terms and have a positive impact on the j values. Among the interactive terms, θ − h₂, θ − h₁, and r − h₂ are almost equally important, whereas r − θ and r − h₁ have a moderate positive impact on j; however, h₁ − h₂ has a negative impact on the j value. The quadric terms are less influential, and all bear a negative impact on j—for instance, θ − θ and h₁ − h₁ bear the lowest significance, whereas r-r and h₂ − h₂ do not contribute much.

For the case of pressure drop, all of the individual terms have a positive impact; the Pareto chart in Figure 7 and p-values in

Table 5. The p-values for f.

Term	p-Value
r	0.0000
θ	0.0000
h1	0.0000
h2	0.0000
r2	0.1232
rθ	0.0000
rh1	0.0002
rh2	0.0000
θ2	0.5765
θh1	0.0000
θh2	0.0000
h12	0.5946
h1h2	0.3063
h22	0.3681

R-squared = 98.7823%. Bold numerals indicate significant factor at 95% confidence level.

illustrate nine significant terms at the 95% confidence level. For instance, h₂ is the most significant parameter, followed by r and θ, having an equal impact and comparatively lesser impact of h₁. The interactive terms r − h₂, r − θ, θ − h₁, θ − h₂, and r − h₁ have a positive impact on f with r − h₂ being the most influential and r-h₁ being the least influential, with h₁ − h₂ as an insignificant term. The quadric terms are non-influential and bear a positive impact except for θ-θ, which has a negative impact.

h₂, r, and θ again have the most effect on pressure drop compared with the leading edge height h₁, in descending order. The obstruction in the flow field is the major cause of pressure drop and these three variables are mainly responsible for this. The leading edge is somewhat aligned with the fluid flow and thus does not play a significant role in the context to the pressure drop. Equations (21) and (22) represent the regression equation for j and f; the insignificant terms are not included in the equation.

$$\begin{gathered} j\left(\times {10}^3\right)=14.1456+0.484726r-0.0771519\theta -1.06379h_1-2.24645h_2-0.213711h_1{h}_2+0.01256r\theta \\ +0.13915r{h}_1+0.299992r{h}_2+0.0200475\theta h_1+0.0207946\theta h_2 \end{gathered}$$

(21)

$$\begin{gathered} f\left(\times {10}^3\right)=407.697-49.9328r-2.7565\theta -41.4562h_1-54.7269h_2+0.257928r\theta +3.12083r{h}_1+4.72761r{h}_2 \\ +0.287108\theta h_1+0.271458\theta h_2 \end{gathered}$$

(22)

Figure 8a shows the 3D response surface plot for j and the effect of θ − h₂ is plotted when the other two variables r and h₁ are at their mean value of 10.5 mm and 2 mm, respectively. The effect of h₂ is more prominent than θ. Figure 8b shows the effect of f—here, again, h₂ has a significant impact on f and θ obeys a linear relationship and has a minimal effect. Figure 8c shows the non-linear variation in j against r-h₁ when the other two parameters θ and h₂ are at their mean values of 52.5° and 2 mm, respectively. Here, h₁ has a more significant effect than r on the heat transfer. Similarly, Figure 8d shows the impact on f—h₁ has a linear relationship, while r has a non-linear relationship. Here, again, h₁ is highly sensitive to pressure losses. To summarize, the impact of θ − h₂ is more significant than r − h₁ on j as well as f.

3.2. Surrogate Modeling with ANN

The ANN consists of neurons imitating the human brain, input layers, output layers (objective function), and hidden layers. It gives weights and biases to the parameters and develops a relationship between the responses without actually digging into the physical mechanism. Thus, it becomes a rapid surrogate model for obtaining the objective function [63]. Figure 9 and Figure 10 show the correlation coefficients (R) of training, testing, and total values as being 0.994, 0.985, and 0.992 for j and 0.999, 0.981, and 0.998 for f, respectively, which are fair figures to indicate a well-trained ANN. In addition, the target-predictor points accumulated near the backward diagonal of the graph, which means that the developed model is quite reliable in predicting the output for different sets of input variables within the test domain. An error histogram is plotted in Figure 11 for j as well as for f training data against the network output; the percentage error lies within a range of ±3% for both j as well as f. The maximum frequency around the zero-error line (highlighted in the red line in Figure 11) ensures an accurate and precise prediction.

3.3. Genetic Algorithm for Multi-Objective Optimization

The genetic algorithm (GA), first introduced by John Holland [64], is a search algorithm based on natural selection and genetics. It has three operators, namely reproduction, crossover, and mutation. Reproduction is based on natural selection based on Darwin’s approach to survival of the fittest. Crossover takes place for the fittest individual in a randomized way followed by mutation, taking into account some of the distinct and important genetic qualities/notions which are again highly randomized. It has the advantage of a random choice for exploitative searches in the design space. The cycle continues until the stopping criteria are reached or the solution is achieved iteratively. The GA is widely used across various research areas; being a stochastic search tool, it allows us to find the global optimal solutions (maxima or minima) based on the objective function. In the present scenario, a multi-objective genetic algorithm (MOGA) was used to maximize the heat transfer and minimize the pressure drop, which are the two objectives of the present study. The MOGA is performed on two surrogate models, viz. RSM and ANN, which use the same DoE. The settings applied for the MOGA are mentioned in

Table 6. Genetic operators and parameters for multi-objective optimization.

Population type	Double-precision floating point vector
Population size	200
Selection operation	Tournament (tournament size equals 2)
Crossover fraction	80% crossover, 20% mutation
Migration operation	20% population moves forward by 20 generations
Pareto front population	35%
Crossover operator	Based on a random weighted average of parents
Stopping criteria	(1) Maximum number of generations = 300 (2) The value of the objective function does not change for the last 100 generations

. The flow chart presented in Figure 12 summarizes all of the steps used in the optimization process. The generated Pareto fronts for RSM-GA (Variant R) and ANN-GA (Variant A) are listed in

Table 7. Pareto front solution points using RSM and multi-objective genetic algorithm.

No.	R (mm)	θ (in Degrees)	h1 (mm)	h2 (mm)	f (×103)	j (×103)
R1	8.109	75	0.712	0.748	48.056	14.443
R2	8.109	74.969	0.718	0.775	48.294	14.492
R3	8.642	74.961	0.781	0.774	49.342	14.908
R4	8.566	74.837	1.014	0.754	51.119	15.112
R5	9.142	74.858	0.751	1.019	53.254	15.595
R6	9.068	74.794	0.871	1.164	55.921	15.95
R7	8.809	74.883	0.987	1.327	57.589	16.164
R8	9.122	74.777	1.062	1.382	61.08	16.578
R9	8.856	74.901	0.757	1.878	62.121	16.839
R10	8.999	74.805	0.92	1.816	64.23	17.036
R11	8.707	74.83	0.975	2.021	64.996	17.134
R12	8.768	74.429	0.814	2.222	66.403	17.343
R13	8.752	74.634	0.733	2.38	67.476	17.494
R14	8.897	74.891	0.77	2.274	67.852	17.519
R15	8.852	74.824	1.077	2.182	69.544	17.614
R16	8.759	74.699	0.902	2.599	72.351	17.976
R17	9.23	74.076	0.881	2.36	73.832	18.088
R18	9.36	74.866	1.02	2.211	74.825	18.119
R19	9.239	73.21	0.935	2.454	75.883	18.278
R20	9.278	72.464	0.865	2.529	76.507	18.356
R21	8.938	74.613	0.783	2.903	77.604	18.537
R22	8.684	74.867	1.092	3.03	79.805	18.637
R23	9.302	74.842	0.83	2.764	80.761	18.788
R24	9.4	74.819	0.866	2.865	84.34	19.095
R25	9.515	74.775	1.13	2.713	86.881	19.205
R26	9.254	74.631	1.042	3.043	87.56	19.341
R27	9.329	74.56	1.109	3.044	89.608	19.484
R28	9.437	74.772	0.985	3.057	89.959	19.548
R29	9.585	74.843	0.918	3.064	91.522	19.683
R30	9.771	74.209	0.945	3.001	93.397	19.792
R31	9.568	74.573	0.985	3.27	96.022	20.042
R32	9.666	74.754	1.174	3.107	97.347	20.054
R33	9.978	74.776	0.782	3.164	98.023	20.189
R34	9.976	74.807	0.829	3.226	100.055	20.343
R35	9.984	74.811	0.969	3.189	101.574	20.409
R36	10.043	74.817	1.18	3.121	104.539	20.539
R37	10.051	74.698	1.071	3.245	105.544	20.667
R38	10.184	74.825	1.056	3.236	107.705	20.803
R39	10.142	74.772	1.25	3.231	109.924	20.901
R40	10.169	74.826	1.365	3.257	113.001	21.071
R41	10.461	74.875	1.249	3.174	115.022	21.179
R42	10.505	74.498	1.262	3.267	118.083	21.389
R43	10.562	74.428	1.267	3.275	119.485	21.47
R44	10.833	74.663	1.163	3.202	121.709	21.559
R45	10.885	74.815	1.277	3.206	125.201	21.738
R46	10.917	74.806	1.29	3.222	126.555	21.815
R47	10.742	74.612	1.534	3.26	128.061	21.88
R48	10.974	74.811	1.299	3.247	128.651	21.936
R49	10.646	74.842	1.812	3.238	130.761	21.946
R50	10.921	74.712	1.473	3.256	130.937	22.036
R51	10.748	74.632	1.742	3.285	132.768	22.089
R52	10.893	74.89	1.735	3.244	135.237	22.202
R53	10.762	74.806	2.036	3.261	138.37	22.279
R54	10.967	74.613	1.824	3.295	139.795	22.428
R55	10.843	74.883	2.18	3.243	142.809	22.445
R56	10.814	74.79	2.212	3.267	143.226	22.461
R57	10.836	74.871	2.3	3.27	145.731	22.553
R58	10.926	74.896	2.272	3.27	147.463	22.652
R59	10.981	74.611	2.383	3.286	151.221	22.788
R60	10.955	74.764	2.474	3.293	152.841	22.832
R61	10.951	74.784	2.558	3.292	154.469	22.87
R62	10.965	74.823	2.609	3.295	156.061	22.923
R63	10.947	74.884	2.749	3.296	158.628	22.974
R64	10.959	74.885	2.86	3.296	161.317	23.038
R65	10.974	74.885	2.86	3.296	161.745	23.059
R66	10.974	74.832	2.996	3.299	164.668	23.114
R67	10.965	74.803	3.086	3.296	166.255	23.129
R68	11	74.839	3.13	3.299	168.325	23.2
R69	10.996	74.96	3.266	3.298	171.378	23.252
R70	10.996	74.96	3.298	3.298	172.069	23.262

and

Table 8. Pareto front solution points using ANN and multi-objective genetic algorithm.

No.	R (mm)	θ (in Degrees)	h1 (mm)	h2 (mm)	f (×103)	j (×103)
A1	8.91	37.963	0.763	0.701	54.032	14.585
A2	9.823	39.472	0.825	0.733	55.068	14.929
A3	9.201	42.812	0.949	0.849	56.34	15.02
A4	9.428	42.577	0.85	1.039	57.8	15.357
A5	9.232	41.039	0.791	1.202	58.394	15.517
A6	9.555	45.229	0.795	1.096	58.939	15.577
A7	9.115	40.535	0.781	1.405	59.682	15.786
A8	9.123	40.536	0.79	1.405	59.728	15.786
A9	9.129	38.568	0.768	1.589	60.649	15.992
A10	9.133	38.569	0.764	1.591	60.663	16
A11	9.225	40.605	0.787	1.631	61.698	16.198
A12	9.14	43.919	0.768	1.666	62.654	16.369
A13	9.327	41.824	0.796	1.847	64.187	16.657
A14	9.169	41.804	0.772	1.963	64.571	16.777
A15	9.305	42.353	0.944	1.986	66.013	16.851
A16	9.234	41.507	0.785	2.273	67.349	17.255
A17	9.382	44.776	0.719	2.176	68.447	17.433
A18	9.39	44.772	0.712	2.176	68.464	17.441
A19	9.41	42.552	0.868	2.386	69.969	17.564
A20	9.462	44.621	0.877	2.378	71.331	17.735
A21	9.433	44.733	0.741	2.536	72.379	18.014
A22	9.525	43.179	0.8	2.697	73.858	18.146
A23	9.441	44.713	0.713	2.73	74.32	18.301
A24	9.236	50.59	0.74	2.66	75.218	18.329
A25	9.311	46.35	0.913	2.919	76.81	18.531
A26	9.311	46.358	0.912	2.919	76.812	18.532
A27	9.46	49.509	0.879	2.769	78.166	18.662
A28	9.509	48.326	0.988	2.84	79.068	18.732
A29	9.52	49.523	1.139	2.847	80.546	18.829
A30	9.287	51.436	0.78	3.163	81.427	19.136
A31	9.582	49.099	0.812	3.075	82.189	19.21
A32	9.638	52.119	0.975	3.181	86.36	19.639
A33	9.39	59.371	1.062	3.262	88.042	19.742
A34	9.787	55.966	0.825	3.225	89.848	20.052
A35	9.888	56.394	1.004	3.206	91.859	20.204
A36	10.008	60.059	0.782	3.238	93.81	20.463
A37	10.036	63.314	0.787	3.281	95.318	20.633
A38	10.145	59.197	0.952	3.263	96.495	20.701
A39	10.316	62.041	0.826	3.266	98.495	20.942
A40	10.302	63.218	1.142	3.273	100.913	21.126
A41	10.479	59.556	1.118	3.285	102.134	21.154
A42	10.416	64.221	1.252	3.224	102.741	21.237
A43	10.633	66.831	0.948	3.27	103.938	21.454
A44	10.337	64.637	1.589	3.274	106.14	21.475
A45	10.816	66.006	0.943	3.273	106.557	21.628
A46	10.74	65.901	1.05	3.296	107.059	21.665
A47	10.782	67.207	1.125	3.288	108.842	21.765
A48	10.928	67.021	1.142	3.283	111.648	21.896
A49	10.821	65.076	1.504	3.282	113.973	21.957
A50	10.663	65.922	1.774	3.282	115.562	22.019
A51	10.856	65.527	1.616	3.285	117.03	22.09
A52	10.757	67.772	1.673	3.287	117.593	22.134
A53	10.878	67.85	1.578	3.288	118.978	22.192
A54	10.652	67.341	1.986	3.283	120.552	22.194
A55	10.966	68.093	1.635	3.297	123.344	22.345
A56	10.831	67.831	1.917	3.283	124.978	22.366
A57	10.859	67.092	2.013	3.285	127.292	22.419
A58	10.959	68.085	1.869	3.285	128.625	22.474
A59	10.964	68.019	1.939	3.298	131.113	22.556
A60	10.935	67.929	2.041	3.297	132.712	22.586
A61	10.872	68.073	2.221	3.282	135.089	22.591
A62	10.975	68.137	2.235	3.295	140.86	22.736
A63	10.999	67.946	2.266	3.297	142.731	22.765
A64	10.966	68.095	2.379	3.294	145.273	22.786
A65	10.981	68.136	2.371	3.295	145.903	22.803
A66	10.972	68.029	2.453	3.298	148.258	22.822
A67	10.98	67.607	2.527	3.3	150.35	22.823
A68	10.993	68.208	2.511	3.295	152.023	22.87
A69	10.991	68.276	2.577	3.295	154.775	22.894
A70	10.999	68.276	2.608	3.3	156.818	22.921

, respectively. The Pareto fronts from both of the surrogate models have a common trait of smaller winglets’ front edge height in comparison to the trailing edge height. The front edge of the winglet is the primary point where the fluid comes into contact with the winglet and it plays a vital role in trapping the flow and redirecting it toward the thermal isolation zone; a larger height would contribute to a higher pressure drop by creating an obstacle in the flow path without much increase in the heat transfer.

3.4. Numerical Investigation of the Selected Models

Figure 13 shows the Pareto fronts obtained from the MOGA using two surrogate models, viz. RSM and ANN. To validate the optimized datasets, a few pareto front points were selected from both datasets—based upon a nearly similar j value—to perform CFD analysis. These points on the Pareto fronts were selected at a regular interval of j values and were further simulated using CFD. Here, Variant A represents the combination of ANN and GA, whereas Variant R represents the optimization based on the RSM and GA. It becomes apparent from the Pareto fronts that Variant A has a higher pressure drop than Variant R until point A21, and beyond that, Variant R has a higher pressure drop. The difference in pressure drop between Variant A and R until point A21 or R16 is quite less than the pressure drop observed beyond these points. Thus, we can say that Variant A has better heat transfer than Variant R for a wider range of design variables. This is because ANN has a better ability to capture the complex relationship between input variables and objective functions.

3.4.1. Flow Field Description

To validate the optimization results, we selected a few pareto front points in a way that the selected points from both variants had nearly similar heat transfer capability. Figure 14 shows velocity streamlines normalized with the inlet velocity (U_in) for the selected models from Variant R and Variant A Pareto fronts and then simulated using CFD at Re = 634.49. Due to constriction in the flow path near the cold tubes, the velocity of the fluid undergoes a substantial increase of approximately three times the inlet velocity. This increase in velocity is caused by flow acceleration, which occurs when the cross-sectional area of the flow path decreases, resulting in a rise in fluid velocity.

The streamlines for the baseline model without winglets follow a straight line without significant mixing (or deflection toward the wake region). As the winglets are attached in the vicinity of the cold tubes, a part of the mainstream flow is diverted toward the recirculation zone. The size of winglets is increased from R6 to R50 and A9 to A50, and with this increase, more of the fluid is redirected toward the recirculation zone. This phenomenon takes place because the winglets’ increasing surface area comes into contact with a large amount of fluid stream and is able to give it a direction.

For a better understanding of the velocity distribution in the flow regime, a normalized velocity contour plot on the z = 0 plane is represented in Figure 15. For the baseline model, an unrestricted high-velocity flow can be observed between the cold tubes and the sidewalls. With the installation of winglets, the flow is directed toward the low-velocity recirculation zones and a better mixing of the fluid stream takes place in the wake region. For the models R50 and A50, the wake region formed is smaller compared to all of the other models due to the larger winglet size.

Figure 16 shows the plots of static temperature normalized with the inlet temperature. The plots show the low-temperature isolation zones formed between the cold tubes for the BLM. Furthermore, for the selected models, the size of these zones decreases with the increase in the winglet size due to the intermixing of the fluid streams. The reduction in temperature is quite evident after the fifth cold tube for models R31, R50, A34, and A50.

Figure 17 shows the turbulence intensity (TI) plots to visualize the velocity fluctuations in the flow regime. It is defined as the ratio of the root mean square (RMS) velocity fluctuations to the mean velocity of the flow—a higher TI refers to higher levels of velocity and fluctuations and it promotes better mixing of fluids. For the BLM, the TI is insignificant, and some fluctuations are observed only in the region downstream of the cold tube walls due to the flow separation taking place past these tubes. In the models selected for analysis, the level of fluctuations increases with an increase in the winglet size and results in a higher heat transfer rate. The wake region is compressed behind the cold tubes for the case of model A9 and is narrower in comparison to R6. For model R16, the spots representing high TI values are formed behind the winglet front edge, while for model A21, high TI spots are formed only near the first cold tube and moderate TI is observed for the successive tube walls. Model R31 experiences higher fluctuations compared to A34, and high TI spots are formed behind the winglets for both of the models—model R31 turbulence is slightly higher behind each cold tube and the wake regions formed are also smaller compared to A34. Similarly, the high TI zones formed for R50 are larger than A50. The higher the TI in the flow regime, the more the disturbance in the flow field, and as a result of this there is better mixing of the hot and cold fluid, and low-temperature zones formed behind the cold tubes are suppressed (cf. Figure 15) by paying the penalty of higher-pressure losses shown in terms of friction factor f in Figure 13.

3.4.2. Performance Analysis of GA and Validation of GA Prediction

Figure 18 shows the percentage deviation in the Pareto fronts of Variant R and Variant A with respect to the simulated values. The deviation in the j values is higher for the design points R6 and R16 compared to A9 and A21, while the deviation in A34 and A50 is on the higher side than in R31 and R50. The maximum deviation is 3.7% for Variant R and 1.7% for Variant A. With the increase in the values of the design parameters, the value of deviation shows a gradually increasing pattern for Variant R while Variant A first shows an increase in the deviation, and then it decreases after attaining the peak value. Now, the deviation in the f values for Variant R is higher for all of the designs than for Variant A.

The maximum deviation is 6.45% for Variant R and 3.12% for Variant A. The deviation in Variant R shows a fluctuating pattern while Variant A first shows a decreasing and then an increasing pattern. Variant A’s overall deviation in j and f is in the acceptable range and shows a better and more precise prediction. This is because the ANN approach has better flexibility and ability to handle complex (non-linear) relationships between input variables and response compared to RSM.

Table 9. Performance comparison of the optimized designs with recent contributions from various researchers.

Researcher		Type of VG	Flow Condition (ReDh)	$\mathit{\eta }$ (%)
Present work (CFD simulations of the selected models)	R6	Curved trapezoidal winglet	634.49	91.62
	R16			99.83
	R31			111.94
	R50			112.28
	A9			95.99
	A21			103.09
	A34			116.40
	A50			117.09
Joardar and Jacobi 1 VG pair [55]		Delta winglet	634.49	30.10
Joardar and Jacobi 3 VG pair [55]		Delta winglet	634.49	36.85
Sharma et al. [65]		Curved VG	200–700	80.0
Naik and Tiwari [66]		Rectangular winglet	2000–4000	10.0
Khan and Li [49]		Rectangular and delta winglet	380–1140	14.4
Gentry and Jacobi [67]		Delta winglet	400–2000	50.0
Sarangi et al. [68]		Curved trapezoidal winglet	300–600	71.0

provides a quantitative comparison between the enhancement factor obtained in the present optimization study for the selected model’s CFD simulation and that obtained from the recent contributions by some of the researchers. The present optimization study offers an enhancement by almost 100% compared to its BLM. The enhancement for Variant A is observed to be greater than for Variant R.

Table 10. Performance comparison of the optimized points with the models used in Joardar and Jacobi [55] with different Re.

Re	HE Variants	CFD		HE Variants	CFD		BLM		1VG		3VG
		j (×103)	f (×103)		j (×103)	f (×103)	j (×103)	f (×103)	j (×103)	f (×103)	j (×103)	f (×103)
280.62	R6	21.17	90.99	A9	21.72	93.69	13.41	79.47	16.69	95.03	17.67	116.07
	R16	23.06	115.29	A21	23.16	108.62
	R31	26.64	134.55	A34	26.24	127.82
	R50	30.35	174.20	A50	29.65	152.00
443.44	R6	18.18	77.77	A9	18.63	79.73	10.19	61.70	12.75	78.33	13.78	96.26
	R16	19.26	93.09	A21	19.68	89.81
	R31	22.80	116.27	A34	22.03	108.38
	R50	25.09	152.19	A50	24.61	134.42
634.49	R6	15.68	59.53	A9	16.20	61.40	8.11	52.17	10.79	72.23	12.59	84.65
	R16	17.60	74.25	A21	17.75	72.53
	R31	20.70	101.21	A34	20.39	90.95
	R50	22.85	135.59	A50	22.38	119.17
850.51	R6	13.00	56.74	A9	13.42	57.98	6.81	44.99	9.87	68.15	11.75	77.89
	R16	14.80	73.77	A21	14.99	70.28
	R31	17.13	93.41	A34	16.39	87.45
	R50	18.02	125.76	A50	17.74	116.10

shows the numerically simulated values of j and f for the selected models and three reference models (BLM, 1VG, 3VG) from Joardar and Jacobi [55] with different Reynolds numbers. The results seem quite interesting. For instance, in the A21 model, with a decrease of 6.7% in the f-value, j increases by more than 42% compared to the 3VG model for the off-design condition. As another example, in the R31 model, with merely an increase of 7.44% in the f-value, there is a significant increase in the j value compared to the 1VG model for the on-design condition. Figure 19 gives a graphical representation of the variation in j values (to the left) and f values (to the right) with different Re. Interestingly, Variant A has f values that are significantly less for all of the Re while the j values are only slightly lesser than for Variant R for models A34 and A50 when compared to models R31 and R50, while model A21 has higher j values and lower f values compared to R16. Model A9 has both j and f values that are higher than model R6. It becomes apparent that Variant A is better than Variant R for on-design conditions as well as off-design conditions. Conclusive results have yielded that the present optimization study holds relevance for both industrial and academic applications. The present study offers a wide spectrum of options to the designers based on their specific requirements and suitability.

4. Conclusions

The present study was undertaken to optimize the performance of the FTHE with curved trapezoidal winglets installed on alternate tubes to maximize j and minimize f. The critical design parameters were first identified before proceeding to DoE using CFD simulation. The radius (r), the angle subtended on cold tubes (θ), the front edge height (h₁), and the trailing edge height (h₂) emerged as the most influential geometric entities for optimization. Considering these parameters, the DoE was prepared using a three-level factorial design which cumulated to a total of 81 design points. The resulting data were used to train RSM and ANN. The R² values for the j and f values were 97% and 90 %, respectively. The deviation in the ANN prediction from the target data was limited to ±3% of the target data. The RSM and ANN data acted as surrogate models for the genetic algorithm (GA) for optimizing the four independent design variables. GA provided two separate Pareto fronts for both surrogate models in a wide prospect.

Some Pareto fronts were selected from the models to evaluate the accuracy of Pareto with respect to the CFD simulations. The maximum deviation in j for Variant R and Variant A was 3.7 and 1.7%, respectively. Similarly, the maximum deviation in f for Variant R and Variant A was 6.45 and 3.12%, respectively. The percentage deviation for Variant A was less than Variant R and hence Variant A is more accurate than Variant R. Additionally, the enhancement factor of Variant A is higher than Variant R.

It was observed that with the increase in the winglet’s size, a larger amount of mainstream flow was diverted toward the thermal isolation zone, enhancing the thermal intermixing, and suppressing the large eddies in that zone. The winglet’s front edge height is always smaller than the trailing edge height.

The CFD study based on the selected models from Pareto fronts confirmed that the heat transfer increased significantly along with the reasonable expense of pressure losses. Thus, the Pareto fronts obtained in this study would be beneficial to academicians as well as researchers across industries. The optimized models were also tested under the off-design fluid flow conditions, wherein their performance was still better for all values of Re.