Natural-convection heat transfer from regularly ribbed vertical surfaces: Homogenization-based simulations towards a correlation for the Nusselt number

Abstract Free-convection heat transfer from vertical surfaces is widely encountered in engineering applications, yet the role played by surface alterations in the heat transfer process and their practical effectiveness are still points of confusion. In this work, buoyancy-driven flows over periodically ribbed vertical plates of different surface micro-textures are investigated, mainly based on an asymptotic homogenization model through which the expensive resolution of the velocity and thermal fields within the inter-rib regions is bypassed, by imposing equivalent effective boundary conditions at a virtual plane surface. Efficiency of the homogenized simulations in detecting macroscopic behavior of the Nusselt number is first assessed, compared with full feature-resolving simulations in which the effects of complex flow patterns, near and within wall corrugations, on the local Nusselt number are captured. Second, the validated model is used to construct a database of numerical results describing deviations of the average Nusselt number over different ribbed surfaces, relative to a corresponding smooth surface. Under the conditions investigated, it is found that surface roughening generally deteriorates heat transfer from vertical surfaces, with slight enhancement for geometries characterized by low thermal slip, for example, rectangular ribs of narrow inter-rib spaces. Finally, a multiple-regression analysis is conducted to formulate a correlation describing effects of the thermal-slip coefficient, the number of ribs, and the Grashof number on the surface-averaged Nusselt number; accuracy of the proposed correlation is attested via further validation. This paper aims to call attention of the heat transfer society to the ability of the homogenization approach to considerably alleviate the computational requirements for relevant simulations and, thus, to significantly accelerate parametric optimization studies.


Introduction and literature review
Buoyancy-induced heat transfer is omnipresent in nature and in technological applications. Being simple in design, cost-effective, highly reliable and free of electromagnetic interference, naturalconvection cooling has gained increasing importance in electronic systems handling low and medium heat loads [1], or even as a pack-up or a supplementary regime in combination with forced convection in the case of high loads [2]. The life and long-term reliability of electronic circuits can be considerably increased by keeping the operating temperature sufficiently below the maximum level allowed by the manufacturer [1,3]. Nevertheless, this has become practically challenging owing to the progressive miniaturization of electronic components, along with the reduction of operational execution time and the increase in power supply [4,5]. Many experimental and numerical investigations have been carried out, seeking enhancement of the heat transfer performance of natural-convection systems through different techniques, including direct immersion into liquids [4,6,7], controlling the spacing of plane heated channels under a variety of boundary conditions [8,9], optimization of heat sink design [10] and, most related to this article, adding surface alterations/extensions to vertical plates exposed to natural convection, which is widely adopted in practical engineering applications such as air solar collectors and cooling of electronics.
In contrast to fins, ribs are by definition relatively small surface disturbances that affect the heat transfer rate from the surface mainly by altering the average heat transfer coefficient, while the increase in the heat transfer area is secondary. Therefore, attaching ribs rather than fins to the heated surface is advantageous in terms of minimizing the overall weight and size of the device, provided that the ribs are proven to be efficient in dissipating heat. This has motivated many experimental and numerical studies investigating the effectiveness of adding ribs, of different shapes and dimensions, to the heated surfaces [2,11,12]. Besides, ribs may be present naturally on the surface, like in the case of electronic components on a chip. Even the detection of adverse effects is important for naturally ribbed surfaces, so that the designer can take this into account while determining the maximum allowable power supply.
When a heated vertical surface is periodically roughened, complex interaction between surface micro-details and the adjacent buoyancy-driven fluid layer takes place, resulting in perturbation of the velocity and thermal fields in the vicinity of the wall, in a quasi-periodic manner, along with macroscopic development of the viscous and thermal boundary layers (cf. [13]). Accordingly, the distribution of the local heat transfer coefficient along the surface deviates both qualitatively and quantitatively from the corresponding distribution along a smooth surface. In particular, the local Nusselt number is significantly affected by mainstream separation and reattachment actions within the inter-rib regions, experiencing successive oscillations along the vertical surface [2,7,11,14]. The overall effect of surface roughening on the heat transfer rate depends on several parameters, including geometric properties of the protrusions, distribution and density of the roughness pattern, thermal properties of the ribs, the Grashof number, the Prandtl number, and the flow regime. However, strict guidelines describing the role of different control parameters are not explicitly available in the literature, a fact that can be deduced from the apparently contradictory findings of previous investigations, with some studies indicating enhancement of the heat transfer performance [12,15], while negative effects of the ribs are reported on many other occasions [11,13,16,17].
In this article, a homogenization-based numerical approach is employed to simulate the macroscale behavior of buoyancy-driven flows over isothermal vertical surfaces, with different periodic micro-patterns of transverse ribs, at a fraction of the cost of full feature-resolving simulations, allowing for generating a large database of numerical results describing the behavior of the average Nusselt number at different geometric and flow conditions, to be analyzed and optimized at a further step in a manner that provides a more mature understanding of complex effects of surface roughening on the heat transfer performance. The present work adopts the effective velocity and temperature boundary conditions constructed and validated by Ahmed et al. [13], extending the scope of previous homogenization studies either not accounting for heat transfer from textured surfaces (cf. for example [18,19]) or limited to low values of the Rayleigh number [20]. In the next section, the physical problem is stated, the governing equations and effective boundary conditions are given, and the geometries under investigation are described in detail. In Section 3, selected cases are considered for simulation with both the feature-resolving approach and the homogenization model, for validation and clarification of the associated phenomena. The generated database of the Nusselt number values for different geometries and densities of patterns is analyzed in Section 4, for optimizing the heat transfer rate from the surface. The role of the thermal-slip coefficient is outlined in Section 5, and a comprehensive equation describing the Nusselt number dependence on the control parameters is formulated and validated. Important notes on the validity and limitations of the adopted physical model are highlighted in Section 6, besides a summary of the main findings of the study.

Outlines of the physical problem
The buoyancy-driven laminar flow over a periodically roughened vertical surface of height (L) and uniform wall temperature at the baseplate and the ribs' surfaces (T w ), higher than the ambient temperature (T 1 ), is considered. Variations in fluid density due to thermal gradient are assumed moderate; the Boussinesq approximation with a linear temperature-density relationship applies (refer to [12,17,21]) so that a temperature-dependent term, representing the buoyancy effect, appears in the momentum conservation equation, allowing for numerical coupling with the energy equation. As a laminar flow is assumed, and sinceT w andT 1 are considered timeindependent, a steady physical model is adopted here; refer to Section 6.1 for further details. In dimensional form, the conservation equations can be written in terms of space coordinatesx i , velocityû i , pressureP, and temperatureT as follows: with q 1 , b 1 ,P 1 the density, the thermal expansion coefficient, and the pressure sufficiently away from the wall, g the gravitational acceleration, and d i1 the Kronecker index (cf. Figure 1 for definition of axes). The dynamic viscosity (l) and the thermal diffusivity (a) are assumed constant; they are defined at the film temperatureT f ¼T w þT 1 2 : Additionally, the Prandtl number (Pr ¼ l q 1 a ¼ a ) and the Grashof number (Gr ¼ gb 1 ðT w ÀT 1 ÞL 3 2 ) are defined, as major parameters controlling natural-convection heat transfer from the surface.

Input variables and optimization objective
Maximization of the average Nusselt number (Nu) over the roughened vertical surface is sought, as the objective of the optimization analysis. The local Nusselt number (Nu) at a point on the surface is defined as: withn the dimensional distance in the profile-normal direction at the point of interest. The conventional definition of the average Nusselt number (Nu conv ) relies on integrating the local Nusselt number over the rough surface and weighting by the total heat transfer area, which is essentially larger than the area of a smooth surface having the same height. Instead, an altered form of the average Nusselt number is defined and used throughout this paper, where the conventional form is corrected by the ribbed-to-smooth surface area ratio, that is Nu ¼ Nu conv Â A ribbed A smooth , to directly indicate the overall deterioration/enhancement of the heat transfer rate in comparison with a corresponding smooth surface, taking into account the increase in surface area due to the presence of ribs; refer to [12]. It can be easily proved that the latter definition of the average Nusselt number (Nu) is equivalent to integrating the local Nusselt number over the rough surface and employing the projected (smooth) surface area as a weight, instead of the ribbed surface area. Accordingly, the average Nusselt number for the case of 2-D ribs can be calculated as: withŝ dimensional distance along the ribbed profile andŜ the total profile length. Four shapes of the roughness element are considered for investigation, specifically, inwardcurved ribs with curves tangent to the base surface, isosceles triangular ribs, circular segments, and rectangular ribs (refer to Figure 1). For all configurations, the geometry of an element, placed in a unit cell of length (l), is characterized by the rib height (e) and width (b). In dimensionless form, the rib height to pitch distance ratio (e=l) and the rib aspect ratio (b=e) are defined. The scope of the present study includes only micro-patterned surfaces, i.e., the pattern periodicity (l) is much smaller than the macroscopic length of the plate (L). This indicates that surfaces have a large number of ribs (N ribs ¼ L=l). The following constraints are set for the geometric variables: For all cases, laminar flow over the plate is considered with values of the Grashof number within the range 10 7 Gr 10 8 , and a constant value of the Prandtl number (Pr ¼ 0.712).

Domain decomposition and homogenized model
Besides the basic feature-resolving numerical approach, that fully captures the surface details and their explicit effects on adjacent flow, the homogenization-based approach constructed and validated by Ahmed et al. [13] is adopted within the present framework. Given the presence of wellseparated length scales, that is, pattern periodicity (l) and plate length (L), the flow field is decomposed into two sub-domains (microscopic/macroscopic) of different characteristic parameters ( Figure 2). Information transfer takes place between the two sub-domains through a virtual matching interface, properly selected for coupling the two problems. The homogenization-based approach eventually seeks a solution to the macroscale problem by imposing upscaled effective velocity and temperature boundary conditions at the matching interface, to be constructed by manipulation of the microscale problem, taking advantage of its amenability to a solution via asymptotic expansions of variables. In brief, the dimensionless conservation equations governing the macroscale problem read: where the macroscopic variables are defined as: using the velocity scale (U ¼ L ffiffiffiffiffi ffi Gr p ). Following the guidelines in [13], effective boundary conditions can be imposed along the vertical plane passing through the crests/outer rims of the ribs (X 2 ¼ 0). The first-order effective conditions, which proved to be efficient for predicting the Nusselt number, read: with ¼ l=L a small parameter characterizing the roughness pattern, R G ¼ ffiffiffiffiffi ffi Gr p the reduced Grashof number, and S 12 the macroscopic shear stress in streamwise direction, defined as: The coefficients (k z , k x , m 12 ) appearing in the effective conditions are all dependent on the rib geometry, that is, rib type, e=l and b/e. They are obtained by solving two-dimensional auxiliary systems of equations through a unit cell of the microscopic sub-domain, in particular, a Laplacelike system for k z and a Stokes-like system for k x and m 12 . The two systems are respectively defined in terms of ad hoc microscopic parameters (h, û i1 , p 1 ) as follows (cf. [13] for derivation): with the boundary conditionsh with the boundary conditions with x i ¼x i =l the microscopic coordinates. Typical numerical solutions of the systems for the case of circular segment (e=l ¼ 0:1663, b=l ¼ 0:4095) are displayed in Figure 3. The coefficients (k z , k x ) are estimated by averaging the fields of (h, û 11 ), respectively, at the plane x 2 ¼ 0, while the coefficient (m 12 ) can be obtained by integrating the field of û 11 over the area below the plane x 2 ¼ 0: The numerical experiments were designed by defining eighty combinations of the geometric parameters (e=l, b=e) for each rib type based on the Optimal-Space-Filling (OSF) algorithm that provides a more uniform extension to the Latin Hypercube Sampling (LHS) method (for more details, the reader is referred to [22]). All numerical simulations were run using Simcenter STAR-CCM þ multiphysics software (versions 15.06.007-R8/16.02.009-R8). The values of the macroscopic coefficients (k z , k x , m 12 ) corresponding to all design points are listed in detail in the Supplementary Material. Response surfaces of the macroscopic coefficients, as surrogate functions of the geometric parameters, are shown in Figure 4 for the different shapes. The Genetic Aggregation algorithm is employed for modeling the response surfaces for all rib types except the rectangular rib, for which the Non-Parametric Regression (NPR) technique is implemented to better handle the high non-linearities in behaviors, as can be realized from the figure. It is clear that, for a given shape, the trends of the three coefficients have some sort of dependence; this will be used in Section 5 to simplify the relationship between these coefficients and the Nusselt number.

Setup of numerical simulations
In this section, three types of numerical simulations, with different targets and approaches, are discussed: (i) simulation of buoyancy-driven flow over a smooth flat plate, to provide a reference for the unperturbed velocity and temperature fields and the value of the Nusselt number relative to which the enhancement/deterioration of heat transfer due to surface roughening can be assessed; (ii) feature-resolving simulations of the flow over ribbed surfaces of different textures, to obtain the detailed perturbed fields and flow structures for a better understanding of the phenomena, and to provide accurate predictions of the Nusselt number for the chosen surfaces; (iii) homogenization-based simulations with first-order effective boundary conditions, to construct the macroscopic fields and the average Nusselt number (Nu) for the different ribbed surface, to be validated by comparing results to the corresponding full feature-resolving simulations. Four cases of flow over ribbed surfaces are considered, that is, the four types of rib geometry shown in Figure 1, with characteristic dimensions e=l ¼ 0:1663 and b=e ¼ 2:4625, are tested. For all cases, the number of ribs is 20 and the Grashof number is 10 8 : The two-dimensional computational domain is sketched in Figure 5, with indication of the imposed boundary conditions. A domain width of 2L is selected, following the specifications by [13], to ensure that the streamwise velocity smoothly vanishes at the far boundary. The uniform pressure values at the inlet and the outlet allow equilibrium with the hydrostatic pressure head in the stagnant flow region. No-slip velocity and uniform temperature boundary conditions are defined along the plane X 2 ¼ 0 for the smooth surface case, while the same conditions are to be imposed on the textured wall when feature-resolving simulations of the flow over ribbed plates are performed. As discussed in Section 2, the homogenization-based simulations require the definition of the effective boundary conditions (Eqs. 6a-c) at the plane X 2 ¼ 0: For the geometric conditions under consideration, the following values of the macroscopic coefficients were numerically estimated: Figure 5. Computational domain and boundary conditions used for numerical simulation. Different categories of boundary conditions are indicated at the wall in accordance with the applied numerical approach, i.e., when ribbed plates are considered, the no-slip velocity and temperature boundary conditions are to be imposed over the real surface structure (feature-resolving approach) or to be replaced by effective boundary conditions at a virtual plane surface (homogenization-based approach). inward À curved Mesh structure and numerical methods are similar to those applied in [13].

Full feature-resolving simulations
The velocity and temperature fields, characterizing the buoyancy-driven flow over the flat plate and the considered ribbed surfaces, are described in Figure 6. The smooth development of the viscous and thermal boundary layers along the X 1 -direction, adjacent to the flat wall, is disturbed by the presence of ribs in a quasi-oscillatory manner. Local increase of the streamwise velocity, accompanied by thermal boundary layer thinning, is realized near the tips/crests of ribs, implying local enhancement of the heat transfer performance. Conversely, the negative streamwise velocity values, clearly captured downstream of the ribs, are attributed to separation eddies where local thickening of the thermal boundary layer is experienced. For a clearer vision, the typical flow patterns near the wall are displayed in Figure 7, revealing complex structures in inter-rib regions where separation and reattachment of the mainstream take place. Proceeding along the plate, the local Grashof number, and the inertia of the flow, increase; consequently, the deflection of the accelerated stream toward the base surface is hindered, reattachment is delayed or even lost, and separation eddies become more pronounced. For the rib height and aspect ratio under consideration, it is noticeable that the recirculation eddies upstream of the ribs are too small to be wellcaptured by the displayed streamlines. Besides the flow patterns, the associated effects on the local behavior of the Nusselt number are presented in the figure. Enhancement of the heat transfer characteristics at tips of the ribs, and deterioration in the recirculation zones are quantified. A secondary peak of the local Nusselt number is detected where the flow approaches or reattaches to the base surface. The Nusselt number macroscopically decreases with the increase in X 1 due to development of the thermal boundary layer (compare levels at low and high local Grashof numbers). The distribution of the local Nusselt number (Nu ¼ @H @X 2 j X 2 ¼0 ) along the flat plate is provided in Figure 8. The numerical results are plotted together with the predictions of Ostrach's similarity solution [23], that is: which perfectly agrees with the present numerical results. Further discussions on the critical assumptions behind the present physical model are presented in Section 6.1, with validation against reference experimental/numerical results from the literature.

Validation of the homogenized model
The homogenization-based numerical model, based on which the optimization analysis will be carried out later in Section 4, is validated here. Inherently, the model targets the macroscopic fields describing the natural-convection flow over the ribbed surface, instead of the detailed patterns of velocity and temperature. Therefore, validation of the model requires comparing its results with the running-average values of the fields resulting from corresponding feature-resolving simulations. For instance, the running-average value of the dimensionless temperature (H) at an arbitrary point (X 1 ¼ a, X 2 ¼ b) is calculated as follows: HðX 1 , bÞ dX 1 : The circular-segment rib is considered for validating the model predictions of the temperature field. As presented in Figure 9 (top), the results of the model are plotted along the virtual interface (X 2 ¼ 0), together with the running-average values of the fully-featured pattern, to explicitly assess accuracy of the temperature effective boundary condition (Eq. 6c). Reasonable agreement in predictions of the thermal slip is realized, i.e., deviation of macroscopic temperature at this virtual plane from its uniform value (H ¼ 1) over an isothermal smooth surface. In addition, the macroscopic temperature profile across a normal section taken at the middle of the plate (X 1 ¼ 0:5) is plotted in Figure 9 (bottom), showing perfect agreement with the running-average profile extracted from the feature-resolving simulation. The profile related to the smooth surface case is also displayed to indicate the deviations experienced when the plate is ribbed.
The main parameter of interest is the average Nusselt number (Nu) over the ribbed plate, corrected with the projected area to directly indicate the total change in heat transfer rate relative to the smooth surface (Eq. 2b). The calculated values of Nu based on feature-resolving simulations of the flow over different ribbed surfaces are given in Table 1. The model predictions are also provided for the purpose of validation, those calculated based on the normal temperature gradient along the matching interface: It is numerically advantageous to exclude the model predictions of the temperature gradient over the first pattern unit (X ¼ 0 ! ), and to replace them in the integration with the corresponding values over a smooth isothermal plate; this can help to avoid numerical uncertainty associated with discretization of the effective temperature boundary condition near the domain inlet, in addition to the model error in this high-convection region (where higher-order effective conditions may be required). It is clear that, for the four cases given in the table, surface roughening negatively affects the heat transfer rate, with the maximum deterioration occurring for inward-curved  ribs. Generally, the model predictions are consistent with the reference values of the featureresolving simulations.

Optimization of heat transfer performance
The rate of heat transfer by natural convection from isothermal vertical surfaces, periodically roughened with transverse ribs of different shapes, characteristic dimensions and pattern densities (refer to Section 2.2 for typical constraints), is the objective of the optimization strategy. In this section, the study is carried out at a constant Grashof number (Gr) of 10 8 ; the effects of variations in the Grashof numbers are discussed later in Section 5. On account of its capability to alleviate the computational requirements of full feature-resolving simulations, the homogenized model which has already been validated in Section 3.2.2 is employed to construct a database of numerical values showing the Nusselt number percentage deviation (D) for different ribbed surfaces (refer to the Supplementary Material). The deviation in the average Nusselt number is measured, taking the smooth surface case as a reference, that is: Again, as discussed at the beginning of Section 2.2, the average Nusselt number over rough surface is corrected (increased) by the ribbed-to-smooth surface area ratio such that (D) directly indicates the total change in the heat transfer rate.
Response surface modeling has been performed to obtain surrogate functions describing the behavior of the Nusselt number percentage deviation (D) with variations in rib characteristic dimensions (e=l, b/e) for different rib shapes and pattern densities ( Figure 10). The same metamodeling algorithms previously chosen for the macroscopic parameters (Section 2.3) are adopted at this step, specifically the Non-Parametric Regression method for the rectangular rib and the Genetic Aggregation technique for other shapes. A comparative look at the levels of D for different rib types reveals that the deterioration in the heat transfer rate decreases by shifting from the blunt inward-curved rib to the more streamlined circular rib (passing by the triangular shape) owing to minimization of surface areas exposed to unfavorable separation eddies (refer to typical flow structures in Figure 7). The relatively favorable effect of increasing the rib aspect ratio (b/e) can be similarly justified. Under specific geometric conditions, the rectangular rib can exclusively promote the heat transfer rate from the surface compared with a smooth plate (D > 0) thanks to its wide outer rim where the main enhancement of the local Nusselt number occurs (refer to the distribution of the local Nusselt number in Figure 7). It is obvious that, for given geometric parameters (e=l, b/e), both deterioration and enhancement of the heat transfer rate become less pronounced by increasing the number of ribs on the vertical surface as the smooth surface case (N ribs ! 1, ! zero) is gradually approached.
Based on the above-mentioned findings, special attention is given to surfaces roughened with twenty rectangular ribs for optimizing the heat transfer rate. For this case, the response surface of the Nusselt number percentage deviation (D) is displayed in Figure 11, as a continuous function of (e=l, b/e). The rib dimensions corresponding to the maximum heat transfer rate were sought; the response surface was analyzed via the Non-Linear Programming by Quadratic Lagrangian (NLPQL) technique, a gradient-based algorithm that implements quasi-Newton methods [22]. A maximum increase in the heat transfer rate of about 1.17% has been attained with the geometric parameters e=l ¼ 0:214 and b=e ¼ 3:741: Overall, for shapes and ranges considered in the present study, ribs either deteriorate the natural-convection heat transfer performance of isothermal vertical surfaces or, at best, enhance it marginally.

Role of thermal-slip coefficient (k z )
In the previous section, the effect of surface roughening on heat transfer performance has been analyzed with variations in rib shape, characteristic dimensions (e=l, b/e), and the number of ribs on the surface (N ribs ). An alternative vision of the previous dependence is elucidated in this  section. Revisiting the conservation equations (Eqs. 4a-c) and the boundary conditions (Eqs. 6a-c) governing the macroscale problem, it can be deduced that the thermal-slip coefficient (k z ) together with the Navier's slip coefficient along the streamwise direction (k x ) and the surface permeability coefficient (m 12 ) can macroscopically mimic the effects of rib geometry on the temperature and velocity fields of the buoyancy-driven flow. Accordingly, the Nusselt number percentage deviation (D) can be expressed as follows: In depth, Eq. (6c) shows that the normal temperature gradient at the virtual interface ( @H @X 2 j X 2 ¼0 ), controlling the Nusselt number, is explicitly related to the thermal-slip coefficient (k z ). However, as the momentum and energy conservation equations (Eqs. 4b and c) are coupled in such a way that the temperature distribution is affected by the velocity field (and vice versa), implicit dependence of the Nusselt number on k x and m 12 is present. It may be advantageous to investigate whether variations in secondary control parameters (k x , m 12 ) can be correlated to the behavior of the main parameter (k z ), so that the function in Eq. (15) may be simplified. The values of the macroscopic parameters for different rib geometries, given in the Supplementary Material, are plotted together in Figure 12; direct relationships between the secondary parameters and the thermal-slip coefficient (k z ) are fairly described via 6 th -order polynomials. Consequently, Eq. (15) reduces to D % f n2 ðk z , N ribs , GrÞ: The dependence of D on the thermal-slip coefficient (k z ) for different numbers of ribs (N ribs ) at a constant Grashof number (Gr ¼ 10 8 ) is presented in Figure 13, based on the numerical results given in the Supplementary Material. Inspection of the figure reveals that, in confirmation to accuracy of the simplified function (Eq. 16), the Nusselt number percentage deviation (D) is well correlated to k z for given values of N ribs and Gr. A critical value of the thermal-slip coefficient, k z Þ critical % 0:025, is indicated; favorable effects of surface roughening on the heat transfer performance of vertical plates are exclusively present below this limit. For both ranges of k z , the magnitude of the percentage deviation (D) decreases when more ribs are present.

Correlation formulation and basic validation
In view of Eq. (16) and the available numerical results, a correlation describing the behavior of the Nusselt number percentage deviation is sought above the critical value of k z , i.e., deterioration in the heat transfer performance due to surface ribbing is considered. The effects of variations in the Grashof number, previously unaccounted for in this article, are analyzed first to be included in the formulation of the correlation at a further step. The values of (D) for vertical surfaces with twenty ribs of different shapes and dimensions were calculated at three values of the Grashof number (Gr ¼ 10 7 , 5 Â 10 7 , 10 8 ) based on macroscopic simulations (refer to the typical results in Supplementary Material); trends are plotted against k z in Figure 14. Notably, the deterioration in heat transfer rate is less considerable at lower values of the Grashof number. This may be attributed to the low level of inertia which allows the mainstream to easily reattach to the base surface so that the unfavorable separation eddies are minimized.
Second, a multiple-regression analysis was performed to formulate a suitable correlation fitting the data given in the Supplementary Material, by implementing the Generalized Reduced Gradient (GRG) non-linear method. Eventually, the following form is proposed: valid at Pr % 0:712 (air), within the following ranges: Predictions of the correlation are plotted in Figure 15 against corresponding numerical values of D, to give an idea about the fitting accuracy.
To check the validity of the introduced correlation for handling rib shapes other than those based on which the correlation has been formulated, the configurations sketched in Figure 16   were considered; the correlation predictions of the Nusselt number percentage deviation (D) for different shapes, at N ribs ¼ 20 and Gr ¼ 10 8 , were compared with the results of macroscopic/feature-resolving simulations, and reasonable agreement was found ( Table 2).
With only the need to calculate k z by solving a two-dimensional Laplace system through a periodic unit cell of the microscopic domain, the present correlation is perfect for preliminary estimation of heat transfer deterioration over ribbed surfaces, before turning to computationally expensive feature-resolving simulations (if necessary). To elaborate this advantage, configuration B was considered for rough estimation of mesh requirements and elapsed wall-clock time for evaluating the Nusselt number percentage deviation (D), via the three methods under discussion. To conduct the full simulation, a total mesh n tot f is used to resolve the features of the whole ribbed surface (with N ribs ¼ 20), and the fields through the full domain (cf. Figure 5), where the conservation equations are discretized. In contrast, the proposed correlation (Eq. 17a) requires to  consider only one representative volume element (cf. Figures 2 and 3), with a grid n tot c , to solve the Laplace system (Eqs. 8a-c) for evaluating k z . The homogenization-based macroscale simulation requires, in addition, to solve the Stokes system (Eqs. 9a-d) for evaluating the coefficients k x and m 12 , so that the effective boundary conditions (Eqs. 6a-c) are available to be imposed on a virtual smooth surface (cf. Figure 5), beyond which the conservation equations (Eqs. 4a-c) are numerically solved, using a grid n tot h : The mesh-dependency analysis in Figure 17 (top) shows that, for a converged solution to be attained, the full simulation requires about one to two orders of magnitude more cells than the homogenization-based simulation and the correlation (to evaluate k z as a prerequisite). Interestingly, the wall-clock time (with six processors) elapsed to solve the Laplace system, so that k z is available to be plugged into the correlation, is five to six orders of magnitude less than the feature-resolving simulation, when mesh-independence is reached; cf. Figure 17 (bottom). Note that the real time required for the homogenization-based simulation (t real h ) includes calculation of the model coefficients and running the macroscale simulation.
6. Comments on physical model validity and concluding remarks

Validity of the physical model
The validity of the physical assumptions adopted in the present work is to be properly assessed, to determine whether the obtained results reflect a true picture of what happens in practice. In  particular, four critical assumptions are to be considered. First, the Boussinesq approximation is reliable under the constraint b 1 ðT w ÀT 1 Þ ( 1 (cf. [12]), such that the density variations due to thermal gradient are significant only in the evaluation of the body force. This assumption, along with the linear temperature-density relationship, is numerically advantageous compared with using the compressible form of the governing equations together with the ideal-Gas law. Second, the dynamic viscosity and the thermal diffusivity of the fluid are considered constant at the film temperature (cf. [14,17]). Third, laminar flow was assumed over the range 10 7 Gr 10 8 (7:12 Â 10 6 Ra 7:12 Â 10 7 in terms of the Rayleigh number). A typical value of 10 9 has been reported/accepted by many investigators as the critical Rayleigh number [2]; however, the criteria based on which the onset of transition is to be determined have been a source of confusion [24]. For instance, Godaux and Gebhart [24] summarized observations from the literature on the value of the Grashof number at the onset of transition (Gr trans: ); accordingly, variations within the range 1:1 Â 10 8 Gr trans: 10 10 (which lies above the limit of the present study, Gr ¼ 10 8 ) can be realized for air natural convection over an isothermal vertical smooth surface. Fourth, the twodimensional flow was assumed; this is reasonable only for the case of transverse ribs (uniform in the spanwise direction), attached to a sufficiently wide plate, over which the flow is laminar. The above-mentioned considerations/limitations necessitate further validation of the present model, mainly against reference experimental results. The classical way for validation relies on plotting the behavior of the local Nusselt number (Nu x ) along the isothermal smooth surface as a function of the local Grashof number (Gr x ), as presented in Figure 18 together with the experimental measurements in [25,26] and the Analytical solution by Ostrach [23]. It should be noted that Nu x is defined based on the vertical position measured from the leading edge, while Nu, plotted in Figure 8, was defined based on the total plate height (L); they are linked as Nu x ¼ Nu X 1 : The local Grashof number (Gr x ) may be calculated as Gr x ¼ Gr X 3 1 , with Gr the plate Grashof number. At relatively low values of the local Grashof number, i.e. Gr x < 1000, the present numerical results significantly deviate from Ostrach's solution, as the boundary-layer assumption used in developing the latter is invalid close to the leading edge of the plate, and applicable only for vertical distances sufficiently larger than the thickness of the boundary layer [23]. Experimental difficulties related to building a perfectly sharp-edged test section may also affect the measurements close to the leading edge [25], resulting in deviations from the numerical values. At relatively large values of Gr x , the present numerical results are identical to Ostrach's solution; however, the experimental results start to deviate, which may be attributed to periodic oscillations close to the trailing edge of the plate and extended effects of room-air turbulence, originally disturbing the outer edge of the boundary layer; refer to [23,27,28].
Providing validation for results of the feature-resolving simulation over a ribbed surface is of much interest. The experimental measurements based on the schlieren technique in [11], with the largest channel aspect ratio, together with the numerical results (using Fluent 6.2) in [17], are used to validate the present fully-featured results of the local heat transfer coefficient (h) over a surface roughened with five square ribs, with the dimensions and thermal conditions shown in Figure 19. With air as the moving fluid, this validation has been conducted at Gr % 3 Â 10 7 , with the fluid properties calculated as discussed in Section 2.1. The results are in reasonable agreement over the edges of the ribs and also through the inter-rib regions. It is notable that a stable thermal boundary layer was optically captured in [11], at the considered value of Gr, which justifies the agreement with the present laminar-flow model.

Conclusions
Natural-convection heat transfer over perturbed surfaces represents a major phenomenon in many engineering and natural processes. Whether the surface is naturally or intentionally roughened, it is vital to analyze the consequent deterioration/enhancement of heat transfer performance for efficiency assessment, design optimization, and determination of the maximum allowable heat load. In this article, laminar buoyancy-driven flows (10 7 Gr 10 8 ) over isothermal vertical surfaces, intensively ribbed (N ribs > 20) with elements of different shapes (inward-curved, triangular, circular, rectangular) and dimensions (height, aspect ratio, etc.) are numerically investigated at a constant Prandtl number (Pr ¼ 0.712). Two numerical approaches, substantially different in scope and computational cost, were tested: (i) feature-resolving simulations with full capability of capturing the micro-structure details of the surface, revealing complex flow patterns including separation and reattachment of the mainstream from/to the baseplate in the inter-rib regions and local thinning and thickening of the thermal boundary layer, besides the associated effects on the local Nusselt number; (ii) homogenization-based simulations with effective boundary conditions, imposed at a virtual vertical plane passing by crests/outer rims of the ribs, to macroscopically mimic the behavior of the flow at a reduced level of mesh requirements and computational cost. For the selected cases, the values of the surface-averaged Nusselt number, based on the two approaches, were compared to validate the model.
Taking advantage of the ability of the validated model to accelerate numerical simulations, it was employed to create a database of results describing the deviation (D) of the Nusselt number, for different rib geometries and values of the Grashof number, relative to the reference value of a smooth isothermal wall. A response surface optimization study was conducted, employing the Non-Linear Programming by Quadratic Lagrangian (NLPQL) technique. It is proved that the presence of ribs generally has a negative effect on heat transfer from the surface, reducing it by more than 10% for some cases, i.e., blunt inward-curved and triangular ribs of large depth. However, a marginal increase in the heat transfer rate has been realized for rectangular ribs of large aspect ratio (width to depth ratio). The maximum heat transfer enhancement, attained under the investigated conditions and ranges, is slightly above 1% and is related to the case of a vertical surface ribbed with 20 rectangular elements characterized by an aspect ratio of 3.741 and a height-to-pitch distance ratio of 0.214, under natural convection with Grashof number equal to 10 8 : It should be highlighted that this marginal improvement is difficult to be validated, as it clearly lies within the experimental uncertainty range (if experimental validation is sought), and even within the numerical uncertainty of the feature-resolving simulations. With the increase in the number of ribs, at a given aspect ratio and fixing the rib height-to-pitch ratio, the positive/negative effects become mild, and the Nusselt number gradually approaches its value over a corresponding smooth surface.
A fundamentally novel view of the contribution of the surface micro-structure to heat transfer by natural convection from vertical surfaces is introduced; the thermal-slip coefficient (k z ), as a single macroscopic control parameter, proved to imitate the inherent dependence of the Nusselt number on surface micro-details under given pattern density and flow conditions. Subsequently, a correlation coupling the effects of the thermal-slip coefficient (k z ), the number of ribs (N ribs ), and the Grashof number (Gr) on the projected-area-based Nusselt number (Nu) was cast based on the available database. The proposed correlation was tested on selected geometries and the results are comparable to those obtained from feature-resolving/homogenized simulations. Accordingly, the correlation is expected to provide a reasonable preliminary estimation of the Nusselt number at an interestingly low computational cost.