Abstract

In the last decade, progresses in additive manufacturing (AM) have paved the way for optimized heat exchangers whose disruptive design will heavily rely on predictive numerical simulations. However, due to typical roughness induced by AM, current wall models used in steady and unsteady 3D Navier–Stokes simulations do not take into account such characteristics. For the development and assessment of novel wall models for AM, a high-fidelity roughness-resolved large-eddy simulation (RRLES) database is built. This article describes the numerical setup and the methodology used for conducting RRLES, from surface generation to postprocessing. In addition, three different cases representing two printing directions plus a streamwise and spanwise isotropic case are investigated. While the roughness distributions are the same in the three cases, the effective slope (ES) is very different, and the impact of this parameter on turbulence and heat transfer is analyzed at different Reynolds numbers.

1 Introduction

In aeronautical engines, heat exchangers have a fundamental role especially in controlling the temperature of the liquid fuel and oil streams. The mass of heat exchangers is also critical to reduce the total weight of the engine. The design of compact and efficient heat exchangers (CHX) is therefore of paramount importance for the future aeronautical engines. Additive manufacturing (AM) offers the opportunity to explore new designs and boost the optimization of CHX [1]. For instance, as shown by Saltzman et al. [2], they obtained an improvement around 10% compared to heat exchangers built with traditional manufacturing.

Nonetheless, surface roughness generated with AM can be substantially larger than with traditional manufacturing techniques. Surface state of metal AM test samples analyzed in several studies have indeed underlined this key point [27]. Moreover, the wall roughness is not isotropic in space and varies according to the main direction of printing [4,5]. Some patterns called welding tracks can appear on the surface [3]. This typical roughness has a significant impact on the performances in terms of pressure drop and heat transfer capacity [2,6,7]. In addition, an important gap between performances obtained from simulations and/or empirical correlations and those measured a posteriori exists.

Many authors have been devoted to the study of the turbulent flow over different roughness types and surfaces experimentally [8] as well as numerically [911]. The presence of surface roughness causes an overall thickening of the boundary layer in addition to the modification of the streamwise structures. Secondary flow motions can also be observed due to spanwise inhomogeneities on the surface, and this induces a better flow mixing. In addition, the transition between laminar to turbulent regimes appears earlier than over smooth surfaces. Concerning the mean velocity profile, a downward shift compared to the log-law profile is observed [12,13]. The latter is commonly called the roughness function (ΔU+). Roughness has also an impact on the Nusselt number since the increase of surface exchange enhances heat transfer [14].

Despite a better understanding of the flow behavior over inhomogeneous rough surfaces, some attempts have been performed to get a universal correlation for evaluating the equivalent sand-grain roughness height ks [15,16,7]. This parameter ks can be expressed as a function of the statistical roughness parameters such as the arithmetic average height (Sa), the root-mean-square height (Sq), the skewness (Sk), and the kurtosis (Ku). The two latter correspond, respectively, to the third and the fourth moments of the height distribution of the surface. Studies have underlined the predominance of the effective slope (ES), defined as the average slope of the surface height along the streamwise direction, in predicting roughness effects [10].
(1)
(2)
(3)
(4)
(5)
where h(x, y) is the surface height, the mean plane is supposed to be at h = 0, and Lx and Ly are, respectively, the streamwise and spanwise lengths.

A well-known threshold equal to ES ≈ 0.35 delimits two regimes: (i) the sparse or waviness regime below this limit and (ii) the roughness or dense one [17]. Nevertheless, such universal correlation that would cover the two regimes has not been found up to now. Furthermore, few articles propose numerical simulation of AM surfaces.

In the aim of building new rough-wall models dedicated to AM surfaces, the creation of a roughness-resolved large-eddy simulation (RRLES) database of representative channel flows emerges as a worthwhile process. The present article describes a methodology to perform parametric RRLES from a limited set of surface parameters to control the surface height distribution and the streamwise and spanwise spatial autocorrelation functions. Section 1 describes the numerical setup, configurations, and the used methodology for performing RRLES. The postprocessing and obtained results are presented and discussed in Sec. 2.

2 Numerical Setup

In this section, chosen configurations and methodology are addressed. In order to ensure periodic channels, recycling conditions are applied on velocity and temperature. This process is particularly exposed here.

2.1 Configurations and Meshes.

First of all, rough surface (RSG) and body-fitted unstructured mesh (RRMG) generators have been developed [18]. The former creates triangulated surfaces stored in the STL format that are used by the latter to compute level set functions. The level set functions are then employed to refine an existing 3D unstructured mesh and to create a conforming body-fitted mesh with controlled cell size, quality, and size gradation. The resulting meshes are suitable for performing RRLES with higher-order finite-volume schemes. The RRMG has been developed within the YALES2 code. YALES2 is an unstructured low-Mach number code developed at CORIA for the direct numerical simulations and large-eddy simulations (LES) in complex geometries [19]. It heavily relies on the parallel volume and surface mesh adaptation.

The RSG/RRMG are utilized to generate rough planar channels but may also be applied to more complex shapes such as cylindrical tubes or plate/fin configurations. RRLES are conducted for periodic channels of size 8H x 3H x 2H in the streamwise (x), spanwise (y), and crosswise (h) directions with H being the half height of the channel, which is equal to 1.0 mm. This domain size was confirmed to be sufficient to have a negligible impact of the periodic recycling on the wall friction [20]. In addition, the half height was selected in order to be close to hydraulic diameters and channel heights encountered for some AM experiments as the L-2x-In case in the study by Stimpson et al. [14].

Concerning the chosen configurations, three different cases representing two printing directions plus a streamwise and spanwise isotropic case have been selected for the analysis (Fig. 1). A smooth channel is also considered as a reference. Roughness parameters that are targeted in this article are the root-mean-square height (Sq), the height distribution skewness (Sk) and kurtosis (Ku), and the ES. While the surface height distributions are the same in the three cases, the effective slope is different due to the different space auto-correlation functions.

Fig. 1
Geometry of considered configurations: (a) perpendicular printing direction to the flow (ES = 0.29), (b) parallel printing direction to the flow (ES = 0.06), and (c) isotropic roughness case (ES = 0.24)
Fig. 1
Geometry of considered configurations: (a) perpendicular printing direction to the flow (ES = 0.29), (b) parallel printing direction to the flow (ES = 0.06), and (c) isotropic roughness case (ES = 0.24)
Close modal

In addition of printing direction, the idea is to generate surfaces that mimic roughness height distribution encountered in AM. An overall range of roughness parameters values can be found in the literature even though there are differences between upward- and downward-facing surfaces [27]. For this purpose, values of the latter for Sq, Sk, and Ku are set, respectively, to 20 μm, around 0.2, and around 4.0. Final values are exposed in Table 1, and probability density functions of height are plotted in Fig. 2 for the three geometries.

Fig. 2
Probability density function of height for the three cases
Fig. 2
Probability density function of height for the three cases
Close modal
Table 1

Roughness parameters of chosen geometries

ESSa (μm)Sq (μm)Sh (μm)SkKu
Front0.2915.620.01870.213.92
Parallel0.0615.620.01810.214.01
Isotropic0.2415.620.02160.213.90
ESSa (μm)Sq (μm)Sh (μm)SkKu
Front0.2915.620.01870.213.92
Parallel0.0615.620.01810.214.01
Isotropic0.2415.620.02160.213.90

With the RRMG, conforming body-fitted meshes are obtained from these geometries. A balance on the element count has to be reached in order to resolve sufficiently the flow while minimizing the cost. In addition, the meshing process has to keep the topology intact. The mesh generation process starts from a Cartesian grid, which is tessellated and adapted numerous times.

The initial numbers of elements in each direction (Nx, Ny, Nh) for the Cartesian grid are fixed in order to initially get isotropic cells of 20 μm. The cell size gradient is equal to 0.1, and the cell size on the rough surface is set to Δyw = 7 μm. The final number of mesh elements for RRLES is about 130 million cells.

All characteristics and performance of meshing process are summarized in Tables 2 and 3. The mesh for the smooth case is a stretched cartesian grid with the dimensionless wall grid spacing of Δx+=13.7, Δy+=5.5, Δh+=[0.5,7.3] in accordance with operating conditions.

Table 2

Meshing characteristics

ParametersValues
Lx, Ly, Lz (mm)8.0; 3.0; 2.0
Initial Nx, Ny, Nz400; 150; 100
Cell size on STL7 μm
Max cell size gradient0.1
ParametersValues
Lx, Ly, Lz (mm)8.0; 3.0; 2.0
Initial Nx, Ny, Nz400; 150; 100
Cell size on STL7 μm
Max cell size gradient0.1
Table 3

Meshing performance

CPURAM/
ESCellsCPUshoursCPU
Front0.29130.1M10244096 h294 Mb
Parallel0.06129.7M10244198 h292 Mb
Isotropic0.24122.9M10245307 h279 Mb
CPURAM/
ESCellsCPUshoursCPU
Front0.29130.1M10244096 h294 Mb
Parallel0.06129.7M10244198 h292 Mb
Isotropic0.24122.9M10245307 h279 Mb

For each configuration, the mesh size is fixed for all the tested Reynolds numbers. The mesh resolution for the present study is discussed in the Appendix. It must be noted that the dimensionless roughness Sayw is only larger than two, which may seem insufficient for the wall resolution. However, with body-fitted grids, this latter has to be sufficient for capturing the wall roughness (curvature and position as illustrated in Fig. 3) and the wall-normal velocity gradient.

Fig. 3
Zoom-in of the mesh on the upper surface for isotropic case (ES = 0.24)
Fig. 3
Zoom-in of the mesh on the upper surface for isotropic case (ES = 0.24)
Close modal

2.2 Statistics and Performance.

For collecting statistics, each RRLES is split into two steps. The initialization step is performed during a given transient time, and then statistics are accumulated. For each step, the durations are summarized in Table 4. For the second step, the dimensionless time t+ defined by t+ = t uτ/H is around 20 in average. Interestingly, due to inhomogeneities of the surface, the time to establish a fully developed turbulent channel flow is reduced compared to smooth cases.

Table 4

Time accumulation for statistics

InitInitStatsStats
Re(ms)No. of FTT (–)(ms)No. of FTT (–)
50000.78.84.050.6
80000.510.72.451.2
17,0000.312.71.251.0
25,0000.0212.70.850.7
InitInitStatsStats
Re(ms)No. of FTT (–)(ms)No. of FTT (–)
50000.78.84.050.6
80000.510.72.451.2
17,0000.312.71.251.0
25,0000.0212.70.850.7

Note: No. of FTT = number of flow-through time.

RRLES CPU costs are presented in Table 5. For clarity, the CPU hours are averaged over all Reynolds numbers cases for each configuration.

Table 5

RRLES performance

ESCellsCPUsCPU hours
Front0.29130.1M56032,800 h
Parallel0.06129.7M56018,500 h
Isotropic0.24122.9M56032,500 h
ESCellsCPUsCPU hours
Front0.29130.1M56032,800 h
Parallel0.06129.7M56018,500 h
Isotropic0.24122.9M56032,500 h

2.3 Numerics and Models

2.3.1 Working Hypotheses.

In each case, incompressible flow simulations are performed. The chosen target bulk Reynolds number range is the fully developed turbulent flow. This is why RRLES are performed at Re = 5000, Re = 8000, Re = 17,000, and Re = 25,000. The fluid kinematic viscosity is set to ν=1.517105m2/s, and the maximal Courant–Friedrichs–Lewy (CFL) number used is equal to 0.8. The wall-adapting local Eddy-viscosity (WALE) subgrid scale model is retained as it is widely used for LES of boundary layer flows [21]. A fourth-order central finite-volume scheme is used, and the four-step fourth-order scheme TFV4A is applied for velocity and scalar transport prediction [22].

2.3.2 Methodology, Boundary Conditions.

The automatic generation of fully periodic channels is very challenging, and ensuring periodic meshes is complex for unstructured grid. The periodicity is though necessary to reach statistical convergence of the flow. Instead of imposing strict periodic boundary conditions in the streamwise and spanwise directions, a Lagrangian recycling method has been developed. This method does not require to add a body force as source term in Navier–Stokes momentum equation to compensate the wall friction. This approach can also be used as a precursor in spatially developing boundary layers by increasing Lout. The idea behind this recycling method is simply to use a time-varying inlet boundary condition with a velocity coming from a distant plane in the channel. This three-step recycling process is based on passive Lagrangian particles as illustrated in Fig. 4. A flowrate is imposed at the inlet and for optimizing performances, and the recycling is done every N time steps, N > 1. A linear interpolation in time is performed between two recycled planes. No rescaling is applied on velocity profile as such profile is unknown a priori. Thus, contrary to what is proposed in the study by Xiao et al. [23], for instance, no target or corrections on velocity field are applied as input of this method. In this article, we select N = 20 as it gives the better compromise between performances and velocity interpolation errors.

Fig. 4
General principle of the developed recycling method
Fig. 4
General principle of the developed recycling method
Close modal

The complete process is described here:

  1. Lagrangian particles are created at the inlet of the grid (spheres at the inlet on the left in Fig. 4) and are translated to the recycling plane.

  2. The target field at the recycling plane, for instance, velocity, is interpolated for translated particles (spheres at the recycling plane in Fig. 4)

  3. This is the relocation step: particles are transferred at the inlet, and the field at this location is updated.

The computational domain is then defined with a recycling zone and a buffer area to avoid any influence of the outlet boundary treatment onto the recycled velocity (Fig. 5). Thus, the location of the recycling plane is primordial and should be set at a given distance from the outlet within the domain. For two parallel planes, this distance was found to be at least equal to the length between both planes, and in our cases, the latter corresponds to 2H. Indeed, perturbations of the velocity field due to the outlet boundary condition may be recycled and injected at the inlet if Lout is below this threshold.

Fig. 5
Computational domain split into two areas: the recycling and the buffer zones
Fig. 5
Computational domain split into two areas: the recycling and the buffer zones
Close modal

In order to decrease interpolation margin errors between the inlet and the recycling plane, both locations should be the equivalent at the wall. Both the generated rough surface and conformal mesh respect this periodicity condition only in terms of wall topology. A no-slip wall boundary condition is applied on rough planes, and the other sides in the spanwise direction are considered as slip walls. This methodology is also applied in the smooth configuration. Validation of the numerical methods may be found in the  Appendix.

Finally, from STL generation to postprocessing, all calculation steps are integrated into a workflow tool. This allows to manage automatically series of runs on a distant super-computer.

3 Simulation Analysis

This section details the different tools used for the subsequent analysis of the RRLES.

3.1 Nondimensional Velocity and Temperature.

For scaling velocity and temperature profiles, a calculated effective distance introduced by Kuwata and Kawaguchi [20] is used. Indeed, due to irregularities of the surface height, this kind of distance is not straightforward to determine as in a smooth wall case. The effective distance is defined in Eq. (6) with hw being the minimal height of the surface. The variable φ corresponds to the xy plane porosity, which is the ratio between the xy plane surface occupied by the fluid and the total xy plane area.
(6)
Computation of the friction velocity uτ (Eq. (7)) is based on the difference between average pressure at the inlet and at the recycling plane as exposed in Fig. 6. The shear Reynolds number is then calculated as Reτ = uτH/ν. and the quantity he+ is defined as he+=heuτ/ν:
(7)
Fig. 6
Principle of the calculation of the friction velocity
Fig. 6
Principle of the calculation of the friction velocity
Close modal

Concerning the friction factor, the Fanning definition f = 2(uτ/Ub)2 is used with Ub the bulk velocity. All these quantities are monitored at each iteration in the LES simulation.

For the analysis of heat transfer, a passive filtered scalar Z¯ is used. This latter can be defined as Z¯=(T¯Tp)/(TTp), where Tp is the temperature imposed at a wall and T is the bulk temperature. This scalar is considered passive, and this hypothesis is valid if the temperature difference has no significant impact on the density, which is assumed here. This is why, the temperature can be replaced by this dimensionless scalar. The equation for this scalar is as follows:
(8)
Thus, Z¯=1 is imposed at the upper wall and Z¯=0 everywhere including the bottom wall. The laminar Prandtl number Pr of this scalar is set to Pr = 0.71, and the turbulent Prandtl number is equal to unity. The diffusivity Dz includes the molecular and turbulent diffusivities.

3.2 Budget Equations.

In the present RRLES, different turbulence budget equations have been computed, dumped, and stored. This is the case for the mean kinetic equation (MKE) (Eq. (9)). Quantities denoted with the bracket 〈〉 symbol are time averaged, and u′ corresponds to the fluctuating component of velocity within the Reynolds decomposition. The quantity νt refers to the turbulent viscosity, which is modeled through the WALE model.

MKE:
(9)
with τijSGS=2νtSij and Sij = (1/2)((∂ui/∂xj) + (∂uj/∂xi)).

4 Results and Discussion

In this section, the wall-unit velocity and temperature profiles as well as the momentum and kinetic energy balance equations described in the previous section are analyzed. This analysis should provide a better understanding of the impact of the wall roughness on the flow and especially the effective slope parameter.

4.1 Impact on Turbulence.

A first qualitative analysis can be done by visualizing the vortices that are generated over the rough surface. The same Q-criterion iso-contours are plotted for the three different cases in Fig. 7 at a Reynolds number of 8000. At a first glance, it is clearly noticed that the parallel case promotes elongated vortices in the flow direction, while the two other cases also feature small vortices trapped between valleys. These trapped vortices are more coherent in the front case due to the orthogonal wavy pattern of the wall. Even if the roughness height distribution is the same, important differences are visually perceptible.

Fig. 7
Q-criterion iso-contour colorized by velocity at Re = 8000: (a) front case (ES = 0.29), (b) parallel case (ES = 0.06), and (c) isotropic case (ES = 0.24)
Fig. 7
Q-criterion iso-contour colorized by velocity at Re = 8000: (a) front case (ES = 0.29), (b) parallel case (ES = 0.06), and (c) isotropic case (ES = 0.24)
Close modal

Mean streamwise velocity profiles are plotted in Fig. 8 for the three cases at Reynolds numbers of 8000 and 17,000 as well as the smooth wall configuration. The impact of ES is clearly visible as the velocity profile is shifted downward with the increasing ES. Differences are also noticeable between isotropic and anisotropic surfaces. Indeed, if the ES in isotropic case is lower than the front configuration, the mean velocity near the wall has a steeper increase for the isotropic case. These results indicate that the orientation of the wall roughness has a nonlinear impact on the mean velocity in the channel. Moreover, increasing the Reynolds number tends to shift downward the profile in the log-law region as expected.

Fig. 8
Velocity profiles (Re = 8000 continuous line, Re = 17,000 dashed line)
Fig. 8
Velocity profiles (Re = 8000 continuous line, Re = 17,000 dashed line)
Close modal

Concerning the fluctuating velocity, streamwise and spanwise components are plotted for Reynolds numbers Re = 8000 and Re = 17,000 in Figs. 9 and 10, respectively. For the streamwise fluctuating velocity urms, and considering the two anisotropic cases, the ES has an impact on the maximum value, and this is latter shifted to the flow stream when changing the roughness direction. This peak is also less pronounced for higher ES value. As for the mean velocity, the isotropic case has a slightly different behavior, and while its ES value is in between the two anisotropic cases, the urms profile is not in between the two of the anisotropic cases confirming the nonlinear behavior of the rough wall with ES.

Fig. 9
Fluctuation of streamwise velocity (Re = 8000 continuous line, Re = 17,000 dashed line)
Fig. 9
Fluctuation of streamwise velocity (Re = 8000 continuous line, Re = 17,000 dashed line)
Close modal
Fig. 10
Fluctuation of spanwise velocity (Re = 8000 continuous line, Re = 17,000 dashed line)
Fig. 10
Fluctuation of spanwise velocity (Re = 8000 continuous line, Re = 17,000 dashed line)
Close modal

For spanwise velocity fluctuations, the maximum of the profile is less affected by the geometry, but its location follows the same trend as the streamwise fluctuations. Increasing the Reynolds number to 17, 000 tends to shift the peaks to the channel center, reduce the sharpness of the peaks for both velocity fluctuations, and slightly increase the maximum value of peaks for spanwise fluctuations.

The reduction of the anisotropy of the velocity fluctuations observed in each case in comparison with the smooth configuration is consistent with previous studies [9,24]. It has been shown that this is caused by the redistribution of the mean roughness wake kinetic energy into turbulence, generating additional normal and spanwise stresses [25]. The streamwise turbulence is also mainly converted into wake energy as the streaks are destroyed by the roughness elements. Increasing the Reynolds number contributes to intensify this phenomenon.

The presence of the wake is slightly visible on the mean velocity profiles when they become negative in the near-wall region for the isotropic and front cases. Mean kinetic energy balance plotted in Fig. 11(compared to Eq. (9) for each budget term) shows that this region is dominated by the work of the pressure drag against the flow direction. Those observations are no more valid for the parallel case since the surface topology does not create significant re-circulation zones, letting the viscous drag to be dominant. Thus, the decreasing of the turbulence anisotropy and the pressure loss are less pronounced for this case.

Fig. 11
Mean kinetic budget terms for Re = 8000 normalized by uτ4/ν from smooth case. 1, ; 2, ; 4, ; 6, ; –7, ; 8, . U0I and U0F refer to mean velocity cross-zero scales for isotropic and front cases, respectively.
Fig. 11
Mean kinetic budget terms for Re = 8000 normalized by uτ4/ν from smooth case. 1, ; 2, ; 4, ; 6, ; –7, ; 8, . U0I and U0F refer to mean velocity cross-zero scales for isotropic and front cases, respectively.
Close modal

For each Reynolds number and all roughness types, the roughness velocity function ΔU+ is studied. The latter is evaluated at he+=100 and compared to the boundary layer log law. Results are summarized in Table 6. The first remark is that ΔU+ is amplified with higher Reynolds number. Values are mainly higher in the front case except at Re = 25,000 where the isotropic case has a larger value. These results point out that the additive manufacturing direction has a strong impact on the mean velocity profile although roughness height distributions are the same.

Table 6

Roughness function values evolution

ReReReRe
ΔU+5000800017,00025,000
Front5.907.8411.113.9
Parallel0.410.691.342.89
Isotropic3.105.4810.215.0
ReReReRe
ΔU+5000800017,00025,000
Front5.907.8411.113.9
Parallel0.410.691.342.89
Isotropic3.105.4810.215.0
Concerning the effect of the effective slope on pressure loss, the friction factor is investigated. The Colebrook correlation [26] (Eq. (10)) is taken as a reference to compare friction factor values obtained in the present RRLES database to the Moody diagram. The definition used in this correlation is Darcy’s one, and the relationship between Darcy and Fanning friction factors is fDarcy = 4f. Therefore, in Fig. 12, friction factors from the Colebrook equation are divided by a factor of 4 for getting comparative results. Two relative roughness Sa/Dh values (Sa/Dh = 10−2 and Sa/Dh = 5.10−3) for Colebrook formula are plotted as the relative roughness of the present cases is intermediate: Sa/Dh = 7.8 10−3.
(10)
Fig. 12
Friction factor values at different Reynolds numbers for each case
Fig. 12
Friction factor values at different Reynolds numbers for each case
Close modal

As expected, for each Reynolds number, the friction factor is higher with the increasing ES. For the front and isotropic cases, values are above the Colebrook correlation ones and tend to slightly increase with the Reynolds number. However, the trend for parallel case is reversed: the curve is below Colebrook expectations and decreases. There is a clear distinction in the friction factor behavior when varying the printing direction. Therefore, this correlation is inadequate for describing the tendency and evaluating a priori the friction factor. It also highlights again the need for better correlations for this type of roughness.

4.2 Impact on Temperature.

The roughness influence on heat transfer is also studied with this RRLES database. Indeed, modifications to the turbulent boundary layer induced by the roughness inhomogeneities can modify the temperature boundary layer and potentially the mean temperature field. As explained earlier, a passive scalar Z¯ is transported and averaged (Eq. (8)). Instantaneous Z¯ fields for the three rough cases are illustrated in Fig. 13. This figure shows that the temperature fluctuations observed close to the wall have wavelengths similar to those of the wall roughness. This is particularly visible in the front and isotropic cases. In the parallel case, the chosen cut cannot represent the wall roughness as the roughness features are aligned with the cut plane. As a result, less short wavelengths are visible in this cut close to the wall.

Fig. 13
Instantaneous temperature field at Re = 8000: (a) front case (ES = 0.29), (b) parallel case (ES = 0.06), and (c) isotropic case (ES = 0.24)
Fig. 13
Instantaneous temperature field at Re = 8000: (a) front case (ES = 0.29), (b) parallel case (ES = 0.06), and (c) isotropic case (ES = 0.24)
Close modal

In order to quantify if temperature transport by turbulence is enhanced by wall roughness, mean temperature profiles are plotted for two different Reynolds numbers in Fig. 14. At Re = 8000, the temperature transport in the boundary layer is enhanced by wall roughness as the mean temperature gradient increases toward a constant value in the whole channel height. The three rough cases have a slightly different behavior close to the wall. It seems that the isotropic case is the configuration that is the closest to the smooth case. For the two anisotropic cases, the presence of coherent structures in the roughness valleys seems to slightly increase the temperature wall gradient and thus the wall heat transfer.

Fig. 14
Dimensionless temperature profiles: (a) Re = 8000 and (b) Re = 17,000
Fig. 14
Dimensionless temperature profiles: (a) Re = 8000 and (b) Re = 17,000
Close modal

For Re = 17,000, the conclusions are different and more intuitive. Indeed, the parallel roughness and smooth cases have almost the same temperature profiles as it would be expected due to the alignment of the roughness elements with the flow. The front and isotropic cases have also the same temperature profile with a clear shift toward the constant gradient imaginary profile. In these cases, the wall roughness promotes temperature transport by turbulence. For this higher Reynolds case, the dependence of the wall heat transfer with the effective slope is better recovered. Again, the different nonlinear behaviors of the temperature with both the effective slope and the Reynolds number illustrate the challenge of finding a Nusselt number correlation for rough walls.

5 Conclusions

This article presents roughness-resolved large-eddy simulations that are representative of the flow obtained in additive-manufactured heat exchangers. The aim of these simulations is to provide a rich database that will ease the development and assessment of rough-wall models. The most challenging task in the building of this database is the generation of representative rough surfaces and conformal unstructured meshes. Three different configurations of parallel rough plane channels with the same roughness height distribution but different effective slopes have been chosen and modeled. The impact of the effective slope parameter, which is directly linked to the alignment of the wall roughness with the flow, has a strong impact on the flow topology, velocity, and temperature profiles, as expected. In these cases, the existing empirical correlations find their limits, and new correlations are needed.

Acknowledgment

The STREAM project has received funding from the Clean Sky 2 Joint Undertaking (JU) (Grant No. 865378). The JU receives support from the European Union’s Horizon 2020 research and innovation program and the Clean Sky 2 JU members other than the Union. This article reflects only the authors’ view, the Joint Undertaking is not responsible for any use that may be made of the information it contains. This work was granted access to the HPC resources of CINES/TGCC/IDRIS under the allocation 2020-A0092B06880 and 2021-A0112B06880 made by GENCI and of CRIANN under the allocation 2012006.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

Nomenclature

f =

fanning friction factor

H =

half height channel

Z¯ =

passive scalar

Reτ =

shear Reynolds number

Sa =

arithmetic average height

Sq =

root-mean-square height

Sh =

maximum distance between peaks and valleys

AM =

additive manufacturing

(C)HX =

(compact) heat exchanger

ES =

effective slope

Ku =

Kurtosis

Pr =

Prandtl number

Re =

Reynolds number

Sk =

skewness

Appendix: Numerical Methods Assessment

We validate our numerical methods by reproducing the classical Reτ = 180 pressure-gradient driven periodic channel flow DNS test case from Ref. [27] (refered to as KMM hereafter). The geometry is Lx = 4πH, Ly = 2πH, Lh = 2H with the same RRLES direction denomination, and the dimensionless mesh resolution is Δx+=8.6, Δy+=8, Δh+=[0.38,5.3]. A first test (T1) is done with the periodic condition methodology, and a second one (T2) uses the recycling boundary condition. For the latter, we set the recycling plane at a distance 2H above the outlet in the longitudinal direction and the input flowrate is chosen to impose the same bulk velocity measured in Ref. [27]. Results are summarized in Fig. 15. Good agreement between tests and references shows that our numerical schemes and recycling boundary condition are appropriate to perform infinite channel flow simulations.

Fig. 15
Velocity mean and fluctuation profiles in the channel flow test case
Fig. 15
Velocity mean and fluctuation profiles in the channel flow test case
Close modal

We discuss now the RRLES mesh size used in this article. Figure 16 compares the results obtained from the Re = 17,000 isotropic configuration with the mesh described in Table 2 (denoted M2) and those obtained with a coarser mesh (M1) whose cell size is 10 μm. It is seen that no significant differences may be observed between the two simulations, on both velocity first- and second-order moments. This suggests that results will be probably weakly improved with a mesh finer than M2, for the analyzed moments.

Fig. 16
Velocity mean and fluctuation profiles for different meshes from RRLES isotropic case, Re = 17,000
Fig. 16
Velocity mean and fluctuation profiles for different meshes from RRLES isotropic case, Re = 17,000
Close modal

Thus, we consider that the dimensionless wall normal resolution y+ obtained with the mesh M2 for this case is a satisfactory reference. Due to the surface alternating of peaks and valleys, making the wall normal mean velocity gradient to vary in space, it is more relevant to observe the statistical distribution of the wall resolution than its mean value alone to evaluate the mesh quality. Figure 17 shows the probability density function of y+ obtained for the considered Reynolds numbers. Obviously, the distributions for the two lowest Reynolds numbers are satisfactory. On the other hand, the distribution for Re = 25,000 is larger than our limit. Nevertheless, it is similar to the one obtained for Re = 17,000 with M1. According to the previous results, we can consider that M2 is adequate for the simulations performed in this article. The analysis for the front and parallel surfaces (not shown here) leads to the same conclusions.

Fig. 17
Dimensionless mean wall resolution PDF from RRLES isotropic cases
Fig. 17
Dimensionless mean wall resolution PDF from RRLES isotropic cases
Close modal

References

1.
Carozza
,
A.
,
2017
, “Heat Exchangers in the Aviation Engineering,”
Heat Exchangers—Advanced Features and Applications
,
S. M.
Sohel Murshed
and
M. M.
Lopes
, eds.,
IntechOpen
,
London, UK
, pp.
149
166
.
2.
Saltzman
,
D.
,
Bichnevicius
,
M.
,
Lynch
,
S.
,
Simpson
,
T.
,
Reutzel
,
E.
,
Dickman
,
C.
, and
Martukanitz
,
R.
,
2018
, “
Design and Evaluation of an Additively Manufactured Aircraft Heat Exchanger
,”
Appl. Therm. Eng.
,
138
, pp.
254
263
.
3.
Cabanettes
,
F.
,
Joubert
,
A.
,
Chardon
,
G.
,
Dumas
,
V.
,
Rech
,
J.
,
Grosjean
,
C.
, and
Dimkovski
,
Z.
,
2018
, “
Topography of as Built Surfaces Generated in Metal Additive Manufacturing: A Multi Scale Analysis From Form to Roughness
,”
Precis. Eng.
,
2018
, pp.
249
265
.
4.
Ventola
,
L.
,
Robotti
,
F.
,
Dialameh
,
M.
,
Calignano
,
F.
,
Manfredi
,
D.
,
Chiavazzo
,
E.
, and
Asinari
,
P.
,
2014
, “
Rough Surfaces With Enhanced Heat Transfer For Electronics Cooling By Direct Metal Laser Sintering
,”
Int. J. Heat Mass Transfer
,
75
, pp.
58
74
.
5.
Simonelli
,
M.
,
Tse
,
Y.
, and
Tuck
,
C.
,
2014
, “
Effect of the Build Orientation on the Mechanical Properties and Fracture Modes of SLM Ti-6Al-4V
,”
Mater. Sci. Eng. A.
,
616
, pp.
1
11
.
6.
Snyder
,
J.
,
Stimpson
,
C.
, and
Thole
,
K.
,
2016
, “
Build Direction Effects on Additively Manufactured Channels
,”
ASME J. Turbomach.
,
138
(
5
), p.
051006
.
7.
Stimpson
,
C.
,
Snyder
,
J.
,
Thole
,
K.
, and
Mongillo
,
D.
,
2017
, “
Scaling Roughness Effects on Pressure Loss and Heat Transfer of Additively Manufactured Channels
,”
ASME J. Turbomach.
,
139
(
2
), p.
021003
.
8.
Jiménez
,
J.
,
2004
, “
Turbulent Flows Over Rough Walls
,”
Annu. Rev. Fluid Mech.
,
36
, pp.
173
196
.
9.
Piomelli
,
U.
,
2019
, “
Recent Advances in the Numerical Simulation of Rough-Wall Boundary Layers
,”
Phys. Chem. Earth.
,
113
, pp.
63
72
.
10.
Napoli
,
E.
,
Armenio
,
V.
, and
De Marchis
,
M.
,
2008
, “
The Effect of the Slope of Irregularly Distributed Roughness Elements on Turbulent Wall-Bounded Flow
,”
J. Fluid Mech.
,
613
, pp.
385
394
.
11.
Busse
,
A.
,
Lutzner
,
M.
, and
Sandham
,
N.
,
2015
, “
Direct Numerical Simulation of Turbulent Flow Over a Rough Surface Based on a Surface Scan
,”
Comput. Fluids
,
116
, pp.
129
147
.
12.
Clauser
,
F.
,
1954
, “
Turbulent Boundary Layers in Adverse Pressure Gradients
,”
J. Aeronaut. Sci.
,
21
(
2
), pp.
91
108
.
13.
Hama
,
F.
,
1954
, “
Boundary Layer Characteristics for Smooth and Rough Surfaces
,”
Trans. Soc. Nav. Arch. Mar. Engrs.
,
62
, pp.
333
358
.
14.
Stimpson
,
C.
,
Snyder
,
J.
,
Thole
,
K.
, and
Mongillo
,
D.
,
2016
, “
Roughness Effects on Flow and Heat Transfer for Additively Manufactured Channels
,”
ASME J. Turbomach.
,
138
(
5
), p.
051008
.
15.
Flack
,
K.
, and
Schultz
,
M.
,
2010
, “
Review of Hydraulic Roughness Scales in the Fully Rough Regime
,”
ASME J. Fluids Eng.
,
132
(
4
), p.
041203
.
16.
Forooghi
,
P.
,
Stroh
,
A.
,
Magagnato
,
F.
,
Jakirlic
,
S.
, and
Frohnapfel
,
B.
,
2017
, “
Toward a Universal Roughness Correlation
,”
ASME J. Fluids Eng.
,
139
(
12
), p.
121201
.
17.
MacDonald
,
M.
,
Chan
,
L.
,
Chung
,
D.
,
Hutchins
,
N.
, and
Ooi
,
A.
,
2016
, “
Turbulent Flow Over Transitionally Rough Surfaces With Varying Roughness Densities
,”
J. Fluid Dyn.
,
804
, pp.
130
161
.
18.
Meynet
,
S.
,
Moureau
,
V.
,
Lartigue
,
G.
, and
Hadjadj
,
A.
,
2021
, “
Automatic Surface and Volume Mesh Generation for Roughness-Resolved LES of Additive-Manufacturing Heat Exchangers
,” ETMM13 Proceedings, ERCOFTAC, Rhodes, Greece, Sept. 15–17, Paper No. 105.
19.
Moureau
,
V.
,
Domingo
,
P.
, and
Vervisch
,
L.
,
2011
, “
Design of a Massively Parallel CFD Code for Complex Geometries
,”
Comptes Rendus Mécanique
,
339
(
2–3
), pp.
141
148
.
20.
Kuwata
,
Y.
, and
Kawaguchi
,
Y.
,
2019
, “
Direct Numerical Simulation of Turbulence Over Systematically Varied Irregular Rough Surfaces
,”
J. Fluid Mech.
,
862
, pp.
781
815
.
21.
Nicoud
,
N.
, and
Ducros
,
F.
,
1999
, “
Subgrid-Scale Stress Modelling Based on the Square of the Velocity Gradient Tensor
,”
Flow, Turbulence and Combustion
,
62
, pp.
183
200
.
22.
Kraushaar
,
M.
,
2011
, “
Application of the Compressible and Low-Mach Number Approaches to Large-Eddy Simulation of Turbulent Flows in Aero-Engines
,” Ph.D thesis,
Institut National Polytechnique de Toulouse
,
Toulouse, France
.
23.
Xiao
,
F.
,
Dianat
,
M.
, and
McGuirk
,
J.
,
2017
, “
An LES Turbulent Inflow Generator Using a Recycling and Rescaling Method
,”
Flow Turbul. Combust.
,
98
, pp.
663
695
.
24.
Shafi
,
H.
, and
Antonia
,
R.
,
1995
, “
Anisotropy of the Reynolds Stresses in a Turbulent Boundary Layer on a Rough Wall
,”
Exp. Fluid
,
18
, pp.
213
215
.
25.
Yuan
,
J.
, and
Piomelli
,
U.
,
2014
, “
Roughness Effects on the Reynolds Stress Budgets in Near-Wall Turbulence
,”
J. Fluid. Mech.
,
760
, p.
R1
.
26.
Colebrook
,
C.
,
1939
, “
Turbulent Flow in Pipes, With Particular Reference to the Transition Region Between the Smooth and Rough Pipe Laws
,”
Inst. Civil Eng.
,
11
(
8
), pp.
133
156
.
27.
Kim
,
J.
,
Moin
,
P.
, and
Moser
,
R.
,
1987
, “
Turbulence Statistics in Fully Developed Channel Flow at Low Reynolds Number
,”
J. Fluid. Mech.
,
177
, pp.
133
166
.