Application of the ‘Nusselt-Equation Simulated Evolution Method’ in Forced Convective Heat Transfer Modeling

Engineers and researchers need proper heat transfer correlations to optimize the design and operation of heat exchangers. These equations are utilized for thermal analysis of the equipment, by applying the Logarithm Mean Temperature Difference (LMTD), Effectiveness–Number of Transfer Units (ε-NTU), Temperature Effectiveness (P-NTU), or any other method. They are a necessity for many project calculations or evaluation of industrial facilities [1-3]. Since the use of wrong or inexact correlations may result in under- or over-sizing of heat exchangers, poor selection of the operational ranges, as well as incorrect estimation of the fouling factors, the use of accurate Nusselt number equations is indispensable to ensure efficiency and profitability of the chemical process industries [4, 5].

Appropriate selection of the analogy among momentum, heat and mass transfer should not be ignored during experimental determination of internal convective heat transfer coefficients in tubular heat exchangers. It not only allows definition of the equation functional form, but also plays a significant role in accurateness of the results. This theory is based on description of the fluid flow through the continuity equations, the Navier-Stokes equations and boundary conditions, hence decoupling the mass species and energy expressions under constant property assumption. The experimental evidences and correspondence amongst both dimensionless expressions indicate that they are similar, which implies that solution of one of the processes conveys to resolution of the analogous phenomenon [6, 7].

Reynolds was the first researcher that realized the link between the energy transport processes. In 1874 he related the heat and momentum transfer for turbulent flows in straight round pipes. Chilton and Colburn improved his postulates in 1934, by extending its applicability to fluids with Prandtl and Schmidt numbers different from the unity $P r \neq 1$ and $S c \neq 1$ [8, 9].

This empirical adjustment to the original approach is often referenced as the Reynolds-Colburn analogy, which resulted in correlations with the mathematical structure of Eq. (1):

$Nu={{k}_{1}}\cdot R{{e}^{{{k}_{2}}}}\cdot P{{r}^{{{k}_{3}}}}$     (1)

where, Nu – Nusselt number; Re – Reynolds number; Pr – Prandtl number; k_i– correlation adjustment coefficients.

Despite monomial power-type equations offers calculation simplicity, they introduce errors in the order of ±25%. In general terms; these correlations are not able to provide truthful approximations over a wide range of experimental data [5, 6, 10].

Prandtl modified the initial model in 1928, by assuming that the boundary layer is comprised by two zones: viscous and turbulent sublayers. He theorized that molecular transport rates prevail on the first region, while a turbulent transport rate does on the second. Correlations derived from this analogy take the Eq. (2) functional form [11]:

$Nu=\frac{\left( f/8 \right)\cdot \left( Re-{{k}_{4}} \right)\cdot Pr}{1+{{k}_{5}}\cdot \sqrt{f/8}\cdot \left( P{{r}^{{{k}_{6}}}}-1 \right)}$     (2)

In Eq. (2) f is the Darcy-Weisbach friction factor. Although Prandtl’s solution offers a significant improvement as compared to Reynolds analogy, because considering the viscous effects, his original equation diverges from experimental data for Pr ≠ 1 or Sc ≠ 1 [6, 7].

A few years later, in 1939, von Karman expanded Prandtl’s analogy by dividing the boundary layer into three zones: viscous, buffer, and turbulent core. He made assumptions alike his predecessor about the relative magnitude between molecular and turbulent transport rates on the viscous and turbulent sublayers, but incorporated the effects over the buffer region by supposing that viscous and turbulent rates were here of the same order of magnitude. The Nusselt number correlation deduced from this analogy is generically represented by Eq. (3) [12]:

$Nu=\frac{{{k}_{1}}\cdot R{{e}^{{{k}_{2}}}}\cdot Pr}{1+{{k}_{5}}\cdot R{{e}^{-0,1}}\cdot \left[ \left( Pr-1 \right)+\ln \left( \frac{5\text{ }Pr+1}{6} \right) \right]}$     (3)

The usage of a three-region model made the turbulent heat and mass transfer theory closer to the realistic flow conditions [13]. However, von Karman’s analogy becomes progressively inaccurate as Pr and Sc increase beyond 30 [6-8].

In spite of their limitations, Reynolds, Colburn, Prandtl and von Karman studies are mandatory references within the context under analysis. In this respect, such analogies are considered the major outcomes of the theory that correlates the heat exchange with friction factor for turbulent convective heat transfer modeling [14].

In most of the methods for characterization of the forced convective heat transfer inside tubes, the correlation functional form needs to be assumed a priori without an exhaustive analysis of the analogy that best describes the experimental data. Even though the mathematical function that relates the Nusselt number with the independent variables is generally selected on the basis of simplicity, compactness and common usage, it cannot be just justified according to these principles. This issue demands special attention, as it can significantly limit the reliability of calculated film coefficients [15, 16].

Despite several parameter estimation methods are applied to obtain heat transfer correlations from experimental data, as briefly summarized by Tam et al. [17, 18], there are only two approaches (according to the authors best knowledge) that eludes anticipated selection of the Nusselt-equation functional form. The first one is symbolic regression, performed through Genetic Programming algorithms, which not only allows determination of the equation constants, but also its mathematical structure [19-21].

However, the resultant expressions are aleatory, have a limited use in practical applications, and are deprived of a theoretical support. The second approach is referred here in as the ‘Nusselt-equation simulated evolution method’, and was devised to link relevant analogies among momentum, heat and mass transfer [18]. It automatically converges into the lowest-errors correlation, through successive modifications of the Nusselt number fitting function. Although physical meaning and generalization of obtained expressions were significantly improved as compared to symbolic regression, this scheme has only been tested versus synthetic data.

Considering previous research gaps, the objective of the actual study was to validate the ‘Nusselt-equation simulated evolution method’ through experimental data, by proposing new forced convective heat transfer correlations and comparing the results with other models that are applicable to single-phase non-laminar fluid flows inside tubes. Key expected contributions are:

• A comprehensive Nusselt number equation, capable of evolving into three different functional forms derived from the Reynolds-Colburn, Prandtl and von Karman analogies;
• Verification of a novel parameter estimation method, applicable to convective heat transfer problems, which does not require anticipated assumption of the correlation functional form;
• A new accurate correlation for transition and turbulent flow conditions, also valid for Prandtl numbers as low as 0.1 (most acknowledged equations are applicable for Pr> 0.5).

The remaining part of this paper is organized as follows: Section 2 explains applied methodology and novelty of the proposed parameters estimation method; Section 3 shows most relevant results, in accordance with the study logical sequence; and Section 4 exposes main concluding remarks.

3.1 Single-range parameters estimation

The full range of experimental data $0.1 \leq P r \leq 10^{3}$ and $3 \times 10^{3} \leq R e \leq 10^{6}$ was used at once, on a first attempt to propose a convective heat transfer correlation for single-phase transition and turbulent fluid flow conditions, inside tubes. Application of the ‘Nusselt-equation simulated evolution method’ provided the following results (see Table 3).

When substituting computed coefficients into Eq. (13), the resultant correlation took the following form:

$Nu'=\frac{0.89\cdot \left( f/8 \right)\cdot \left( Re-136.2 \right)\cdot Pr}{1+10.478\cdot \sqrt{f/8}\quad\cdot \left( P{{r}^{2/3}}-1 \right)}$     (20)

In spite of the strong linear association that was determined through Eq. (20), by comparing the Nusselt number calculated values (Nu') against the reference ones (Nu), the more accurate results were reported for Eq. (6) suggested by Taler. In contrast, power-type correlations like Taler’s Eq. (7) were confirmed less appropriate to reproduce experimental results over a broad range of the independent variables (see Table 4).

Table 3. Single-range parameter estimation results

Notes: ¹ Stop criteria:

I. Maximum number of generations exceeded;

II. Average change in the fitness value less than options.

where, f_obj – objective function value once the algorithm stops; i – number of iterations.

Table 4. Assessment of Eq. (20) vs. benchmark correlations

where, R² – coefficient of determination; e_ave – mean relative error,%; e _max – maximum relative error,%; SSE – sum of squared errors. Utilized statistical and error metrics are described in Appendix A.

A 3D plot of the relative errors introduced by Eq. (20), as function of Prandtl and Reynolds numbers, shows larger deviations for $\operatorname{Pr} \leq 3$ and $R e<10^{4}$, with the first variable exerting the greatest influence (Figure 2).

3.2 Two Re-intervals parameters estimation

Regression by intervals is a common practice during experimental determination of convective heat transfer coefficients. In this respect, authors as Camaraza-Medina et al. [31, 32] performed parameters estimation for different ranges of the Reynolds number, comprehending that the flow regime have a direct impact on heat transfer intensification. Besides, some of the theories that correlate momentum, heat and mass transfer for turbulent flows are based on a velocity profile approach [34-36].

2.png

Figure 2. Relative errors 3D plot for Eq. (20)

Taking the above into consideration, a second approximation of the experimental data was carried out by batches, selecting two intervals of the Reynolds number: $3 \times 10^{3} \leq \operatorname{Re}<10^{4} \quad$ and $\quad 10^{4} \leq \operatorname{Re} \leq 10^{6}$. Parameters estimation results are shown below (see Tables 5 and 6).

Substitution of computed coefficients into Eq. (13) provided the following expression:

$Nu'=\frac{{{y}_{1}}\cdot \left( f/8 \right)\cdot \left( Re-{{y}_{2}} \right)\cdot Pr}{1+{{y}_{3}}\cdot \sqrt{f/8}\quad\cdot \left( P{{r}^{2/3}}-1 \right)}$     (21)

where,

$[{{y}_{1}},{{y}_{2}},{{y}_{3}}]=\left\{ \begin{align}  & [0.9052,\text{ 7}\text{.80},\text{ 1}0.752],\text{ if }3\cdot {{10}^{3}}\le Re<{{10}^{4}} \\ & [0.8902,\text{ 243}\text{.75},\text{ 10}\text{.478}],\text{ if 1}{{\text{0}}^{\text{4}}}\le Re\le {{10}^{6}} \\\end{align} \right.$

Table 5. First Re-interval parameters estimation results

Even when Eq. (21) provided slightly better results as compared to Eq. (20), accuracy of Taler’s Eq. (6) remains the higher to this point (see Table 7).

Table 6. Second Re-interval parameters estimation results

Table 7. Assessment of Eq. (21) vs. benchmark correlations

3.png

Figure 3. Relative errors 3D plot for Eq. (21)

3D plot of the relative errors introduced by Eq. (21) still exposes sudden deviations for low Prandtl numbers (see Figure 3). Potential reasons behind this outcome are the influence of this variable over the exponents of the Reynolds and Prandtl numbers, as well as the distinctive thermal behavior of liquid metals and analogous fluids (Pr<< 1).

As confirmed by Gnielinski [37] and Taler [10], power exponents at the Reynolds and Prandtl numbers strongly depend on the Prandtl number in turbulent flow. If Pr is smaller than one, power exponents with Re and Pr are close to each other. Besides, Taler highlighted that the Re exponent increases with the Prandtl number, while the Pr exponent responds inversely.

On the other hand, in low Prandtl-number fluids the contribution of the molecular thermal conduction to the total heat transfer is higher than other streams. Since having high thermal conductivity coefficients, thus very low molecular Pr that result in a much thicker thermal boundary layer as compared to the hydrodynamic one, it becomes difficult to correlate them concurrently with moderate and high Prandtl-number fluids [6, 38].

3.3 Two Pr-number intervals parameters estimation

Despite a few researchers tackled heat transfer problems by examining different ranges of the Reynolds number, others like Dalkilic & Wongwises [39], Gnielinski [23] and Taler [10] proposed film coefficient correlations based on Prandtl-number intervals. Likewise, it is noteworthy that accurate models like Polley’s [40], Sleicer-Rouse’s [41] and Sandall et al. [30] are supported on correlations where the Prandtl number have a significant influence as independent variable.

Therefore, and according to this study methodology, next approximation of the experimental data was carried out by batches over the Prandtl number: $0.1 \leq \operatorname{Pr} \leq 3$ and $3<\operatorname{Pr} \leq 10^{3}$. Parameters estimation results are summarized below (see Tables 8 and 9).

Table 8. First Pr-interval parameters estimation results

Table 9. Second Pr-interval parameters estimation results

Substitution of the above coefficients into Eq. (13) produced the following Nusselt number expression:

$Nu'=\frac{{{z}_{1}}\cdot \left( f/8 \right)\cdot \left( Re-{{z}_{2}} \right)\cdot Pr}{1+{{z}_{3}}\cdot \sqrt{f/8}\quad\cdot \left( P{{r}^{2/3}}-1 \right)}$     (22)

where,

$[{{z}_{1}},{{z}_{2}},{{z}_{3}}]=\left\{ \begin{align}  & [0.9713,\text{ 205}\text{.05},\text{ 12}\text{.9}52],\text{ if }0.1\le Pr\le 3 \\ & [0.8761,\text{ 147}\text{.30},\text{ 10}\text{.300}],\text{ if 3}<Pr\le 1000 \\\end{align} \right.$

The resultant convective heat transfer correlation, given by Eq. (22), provided a better fit to the reference data (Nu) as compared to benchmark expressions (Table 10). Unlike the other equations, relative errors were lower than 20% for all data points (Figure 4).

4.png

Figure 4. Relative errors 3D plot for Eq. (22)

Table 10. Assessment of Eq. (22) vs. benchmark correlations

The constant c₄ in Eq. (13) is related to the dimensionless velocity at the hypothetical distance from the tube surface to the boundary limit between the viscous (laminar) sublayer and the turbulent core. It depends on the viscous sublayer thickness and was expected from c₄ = 5 for the three-sublayers model of von Karman to c₄ = 11.7 for the two-zone boundary layer theory of Prandtl [22]. Despite slightly outside of the theoretical range, the results from this study c₄ = z₃ ={12.952; 10.300} were close to the coefficient 12.475 found by Taler [25]. In both cases, the additive form of the denominator arose from the two resistances in series of the Prandtl model for momentum and energy transfer [35].

It was remarkable that all runs of the ‘Nusselt-equation simulated evolution method’converged into a correlation that corresponds to Prandtl analogy. This emphasize that monomial power-type correlations based on Reynolds analogy do not provide good approximations, since assumed that the turbulent transport rates were much greater than the molecular transport rates, neglecting the later. The viscous effects are significant and can make comparable or even greater contributions than the turbulent effects in the near-wall region. It also underlines that Nusselt-number functional forms according to von Karman analogy are neither accurate if the studied data points include high Prandlt numbers. Deviations beyond Pr= 10 were attributed to the assumption of a fixed boundary between this model sublayers, independent of Reynolds and Prandtl numbers [6-8, 36].

3.4 Comparison with other correlations

Considering previous results, Eq. (22) is recommended for calculation of the mean convective heat transfer coefficient for single-phase transition and turbulent fluid flow conditions, inside tubes. Besides accurately describing the experimental database used for this study, computed Nusselt number values agreed very well with those obtained through six known correlations, those of Skupinski et al, Seban & Shimazaki, Gnielinski, Camaraza-Medina et al., Nalavade et al. and Sandall et al (see Figure 5). The selected models are the widely-acknowledged correlations in the technical literature [26-32].

5.png

Figure 5. Scatter plot for computed Nu values

The relative deviations amongst Nusselt number values, calculated by Eq. (23), were used as a performance metric for quantitative comparison of the suggested correlation versus the reference ones:

${{d}_{Nu}}=\frac{\left| Nu''-Nu' \right|}{Nu''}\text{ }\cdot 100\text{ }%$     (23)

where, d_Nu – relative deviation over the dependent variable, %; Nu'– Nusselt number value calculated by means of Eq. (22); Nu'' – Nusselt number value calculated with the selected models.

Owing that mean divergences were lower than 12% (see Figure 6), besides not exceeding 20% differences in 92.59% of the data points (see Table 11), this research correlation is considered appropriate for academic purposes and practical applications within the studied validity ranges: $3 \times 10^{3} \leq R e \leq 10^{6}$ and $0.1 \leq \operatorname{Pr} \leq 10^{3}$. 

6.png

Figure 6. Relative deviations amongst Nu values

Table 11. Relative deviations percentage distribution

The largest deviations were found when comparing Eq. (22) versus the results from Gnielinski equation, since heat transfer modeling of the transitional flow zone is not as predictable as the turbulent regime. At one side, the lack of experimental data is most likely the reason for limited understanding and little design information about the transitional flow region, so it is yet considered a metastable and complicated zone [42]. At the other side, contemporary review of Gnielinski equation indicated a considerable amount of scatter (at least ±20%), presumably due to the antiquity of the data he used and the absence of modern instrumentation [43, 44].

Eq. (22) can be generalized to a tube of finite length and streams having large property variations due to remarkable temperature differences, by introducing two well-known correction factors [25, 37]:

$Nu'=\frac{{{z}_{1}}\cdot \left( f/8 \right)\cdot \left( Re-{{z}_{2}} \right)\cdot Pr}{1+{{z}_{3}}\cdot \sqrt{f/8}\quad\cdot \left( P{{r}^{2/3}}-1 \right)}\cdot \left[ 1+{{\left( \frac{d}{L} \right)}^{2/3}} \right]\cdot {{\left( \frac{Pr}{P{{r}_{w}}} \right)}^{0.11}}$     (24)

A correction based on temperature rather than the Prandtl number ratio (right-end multiplier) is recommended for gases [25, 37]:

${{K}_{\mu }}={{\left( \frac{{{T}_{b}}}{{{T}_{w}}} \right)}^{a}}$     (25)

where, Pr_w – Prandtl number evaluated at wall temperature, T_b – bulk temperature, K; T_w – wall temperature, K; a – correction factor exponent, equals to 0 for cooling and 0.45 for heating. While the last value was correlated by Gnielinski for 0.5 <(T_b/T_w)<1.0, there are other recommendations in the literature [45].

Despite the boundary condition assumed to gather the experimental database, recommended correlation can be used for either uniform wall heat flux or uniform wall temperature applications, if Re>10⁴. While in the laminar and transitional flow regimes the boundary condition has a significant influence on the heat transfer characteristics, turbulent flow is practically insensitive to the different conditions [46-52].