Research on Aerodynamic Roughness[1]
Mu Qingsong1,2 Miao Tiande1 Dong Zhibao2
(1Department of Mechanics, School of Civil engineering and Mechanics,Lanzhou University, Lanzhou, 730000, China; 2Key Laboratory of Desert & Desertification , Cold and Arid Regions Environmental and Engineering ResearchInstitute , Chinese Academy of Sciences , Lanzhou 730000 , China)
Abstract Theoreticalanalysis and numerical simulation are used to explain important experimental phenomena.For examples, when thewind velocity above the boundary layerincreases,the aerodynamic roughnesslength decreases nonlinearlywith increasing wind velocity, but shear velocity increases proportionally tothe increase wind velocity, and when the distribution density of roughness elements is artificially increased, both roughness length and shear velocity increase at a progressively slower rate. The theoretical resultsclosely paralleled the experimental data,and contradicted the traditional opinion that the roughnessonly relates to thegeometrical shape of a rough wall.Recomposed Nikuradse’s traditional hydro-dynamical theory of roughness, which does not consider the influence of the distribution density of roughness elements on roughness length and shear velocity,our theoretical resultssuggest that the dynamical natures of roughness length and of shear velocity aredetermined by both thegeometrical shape of the rough wall and the fluid field around the wall.
Keywords:wind velocity; roughness length; shear velocity; displacement height; geometric shape of roughness elements; Nikuradse’s hydro-dynamical theory of roughness
1. Introduction
Turbulent flows in the boundary layers, such as wind blowing over soils, deserts, shrubs or gravel fields, will cause erosion. Roughness elements such as clods of soil, stones and plants at the surface can partially protect that surface from erosion, and thus offer important protection for the environment. In a turbulent boundary layer, the essential dynamical parameters are the shear velocity and aerodynamic roughness length; as result, the relationship between these parameters and wind velocity above the boundary layer, and the relationship between these parameters and the shape, volume, and distributing density of roughness elements, have attracted the interests of geologists and physicists. For nearly seven decades, researchers have been reporting the results of fieldobservations andlabexperiments including naturalwalls(for instances, nude-soil, sand-bed and gravel-surface) and artificial rough walls(for instances, array of hemi-spheres, array of columns and array of cones). The main results of these researches can be summarized as follows: under a certain wind velocity, the turbulent field near a rough wall is respectivelydefined isolated-roughness flow, wake interference flow and skimming flow,which respectivelycorresponds to sparse, dense and more dense distributing density for the roughness elements.
Some researchers, such as Morris[1], Lee and Soliman[2], Wolfe and Nickling[3-4], consideredthe roughness length of a rough wall as a geometrical parameterthat depended only on the volume, shape and distributing density of its roughness elements, and concluded based on their experiments that roughness length increases nonlinearly with increasingdensity; however, the opposite occurs when the density exceeds a specific threshold value. Others, including Donget al.[5], considered roughness length to be an aerodynamic parameter that is relevant to the flow field near the wall besides the geometrical shape of the wall. Many studies ofshear velocity have concluded that the relationship between shear velocity and the distributing density of roughness elements is similar to that between roughness length and this density, and that the shear velocity increases nonlinearlyas the wind velocity above the wall regionincreases and progressively becomes proportional to the wind velocity. However, Bagnold[6] promoted the idea that shear velocitywas only determined by the flow field above the turbulent boundary layer.
In our review of the experimental and theoretical studies that have been performed in this research field, we found that different scholars focused on different aspects of the dynamics of roughness length and shear velocity,and this difference in focus revealed some importanttheoretical problems that have not yet been solved.These include the following questions: Is roughness length a geometrical parameter that is only determined by wall shape or is it an aerodynamic parameter related to both wall shape and turbulent flow-field near the wall? Is shear velocitydetermined only by the flowfield above the wall region, and does it have no correlation with the geometrical shape of a rough wall? In this paper, based on our previousresearch[7],we focus on these problems and develop an improved theory that clarifies how the density of arbitrary roughness elements and wind velocity influence roughness length and shear velocity, and present some new viewpoints about the dynamical characteristics of roughnesslength and shear velocity.
2. Logarithmic profiles of wind velocity around rough walls
We begin our analysis by assuming an ideal rough wall that the roughness elements of uniform size and shapedistribute on a horizontal solid surface (Fig.1), and we takea coordinate system withthe origin on the surface, zaxis vertical above the surface and xaxis in the wind direction. According to Prandtl's mixing-length theory, the momentum equation of air in a steady, one-dimensional, incompressible flow can be expressed as follows
， (1)
Where *is the shear stress on the wall due to wind, zxis the shear stress, gis the atmospheric density，lis the mixing length，and u isthe average velocity of the turbulent wind.In the turbulent boundary layer near a smooth wall，there is a thin laminar layer that clings to the solid surface, and above this layer lie a transitional region and a turbulent region. Prandtl[2] considered that the wind velocity (u) in the turbulent region should satisfy the momentum equation (equation 1), and proposed that the mixing lengthl=z, where =0.4is the karman constant.The center of wind’s shear stress (*)on the smooth wall occurs at z=0 because* is only caused by the wind’s friction force on the wall.However, the center of*occurs atz=d on the roughness wall (where d is the displacement height) because * consists of the drag force exerted by roughness elements on the wind and the wind’s friction force on the uncovered part of the wall.Thom[8] and Jackson[9] proposed that the mixing length(l)should be assumed to be directly proportional tothe difference between the actual height and the displacement height (i.e., z-d) in the turbulentlayer near a rough wall instead of being proportional tozin the turbulent layer near a smooth wall. The turbulent velocity can be obtained by substituting (z-d)for l in Equation (1),then integrating the resulting equation with respect to z, giving the following equation
（2）
Fig. 1 Schematic diagram of the wind velocity profile in a turbulent boundary layer.
In which u*=(*/g)1/2is shear velocity, and z0 is an integral constant, named the aerodynamic roughness length. We use to expressthe layer thickness, when z ≥ ,assuming that the wind velocity becomes uniform and equals to U.. From this assumption, we can derive
（3）
3. The shear velocity of a rough wall
Based on dimensional analysis[8-10], the drag force (F) exerted by one roughness element on the wind can be expressed as follows
（4）
Where h is the vertical height of roughness elements, is the horizontal width of the roughness element in the direction of wind velocity, Cris the drag coefficient,and Ais the coefficient of the upstream projected area. The first part that contributes to the wall’s shear stress*is the stressrcaused by drag forces exerted by the roughness elements, and can be written as
（5）
Where n is the number of roughness elements per unit area of the wall(i.e., the density of the roughness elements). The second part that contributes to the wall’s shear stress*is the stresss caused by the wind’s friction force on the uncovered part of the wall. Thus, we have
（6）
Whereis the area fraction of uncovered part on one unit area of the wall.scan be expressed as
（7）
Where Cs is the drag coefficient of the uncovered part of the wall. We define the roughness concentration()as the sum of the upwind projected areas of the roughness elements to the wall area
（8）
and define the Reynolds number (Re)of the roughness elements as
（9）
In which is the kinematical viscosity coefficient of air. In nature, Cr and Csare a function of Re and Many field observations and laboratory experiments (e.g., Marshall[11], Lyles and Allison[12], Gillette and Stockton[13], Musick and Gillette[14]) have concluded that Crand Cswill decrease with increasing Rebut will become approximately constant for very large values of Re. in addition, with increasing , Cr will decrease due to the increased interference between roughness elements and Cs will decrease due to the increased area of the wall that is covered by roughness elements.
Suppose that Cr0 is the drag coefficient of the roughness elements in theisolated-roughness flow, and we can ignore the interference between roughness elements because of their sparse density. In this case, we can say that Cr0 depends only on Re. Whentends to zero, we assume that Cris a linear function of, then
（10）
In which r is a proportionality coefficient. Whenbecomes larger, we divideinto m parts, each of which ()equals/m. If the roughness concentration is only , then the drag coefficient Cr1 can be express as
（11）
Let us assume that the relationship between Cr2 and Cr1 is the same as the relationship between Cr1 and Cr0 when the roughness concentration is2, namely
（12）
From equations(11) and (12), we have
（13）
When the roughness concentration is m, Cr equals
（14）
As mapproaches infinity, equation (14) becomes
（15）
If we assume that the same analysis can be performed for Cs, we obtain
（16）
Wheresis a proportionality coefficient, and Cs0is the drag coefficient for a smooth wall with no roughness elements. As Re approaches infinity, Cr0 and Cs0 becomeconstants. This produces the following result
（17）
From equation（17）, we know thatr increases with increasing for small values of but decreases with increasing for large values of , butsalways decreases with increasing . The experimental study of Schlichting[10]and Marshall[11]concluded thats becomes much smaller thanr whenexceeds a threshold value, which Schlichting and Marshall determined lies within the interval [0.003, 0.1].
4. The dynamical nature of roughness
If we ignore s in*, according to equation（3）and equation (17), we have
（18）
A result of d=0.7h was obtained in experimental studies by Blihcoand Partheniades[15]and Kamphius[16]. Because decreases with increasingReandRe=Uh/v anddecreases with the decreasing, we can conclude from equation (18) that will decrease when increases, which is in accordance with experiment data[5].
If ∂z0/∂=0, then we have=r-1, a critical number. For r-1, z0 increases withincreasing but forr-1, z0 decreases withincreasing . As a result,z0 has its maximumat = r-1.
Suppose that the geometrical shape of the roughness elements is certain, the changeinis only caused by the changein .In this case,z0 increases with increasing n when nA-1h-1-1r-1,but decreases with increasing nfor nA-1h-1-1r-1, in agreement with the experiments of Morris[1], Lee and Soliman[2], and Wolfe and Nickling[3-4].
5. The dynamical nature of shear velocity
Suppose thatscan be ignored. Equation (17) then becomes
（19）
In the hydro-rough situation, Cr0remains approximately constant,so u* is proportional toU when Reincreases. If we let ∂u*/∂=0,we also define the critical point = r-1,at which u*reaches its maximum.The relationship between u* and is similar to that between z0 and,and consequently, the relationship between u* and n is also similar to that between z0 and n.When is close to zero, u* increases rapidly with increasing ,but when becomes larger,u* gradually begins to increase more slowly with increasing and finally begins to decrease with increasing .The relationship between u* andfollows a typical course in which the isolatedroughness flow gradually becomes the wake interference flow, and the wake interference flow finally becomes the skimming flow as increases within the turbulent boundary layer, corresponding to the results of experiments by Marshall[11], Lyles and Allison[12], Gillette and Stockton[13], and Musick and Gillette[14].
Experiments by some researchers[5] have found that u* and z0 always increase with increasing and finally tend towards a constant value for large values of ,but have not reported that u* and z0 decrease with increasing after exceeds a certain critical threshold.One possible explanation is that the critical numberr-1 in these experiments was larger, so that could not exceedr-1 even if the cover coefficient, defined by the sum of the vertical projected areas of the roughness elements to the wall area, approaches one.
5. Comparing theoretical models of the dynamical nature of roughness length and shear velocity with experimental data
We used the data obtained by Dong et al.[5]to test the theoretical equations developed in previous sections of this paper.In their experiment, thin cylindrical sticks with a height of 16cm and a diameter of 0.3cm were chosen to create arrays of roughness elements on a 0.52.0 m plate attached to the bottom of the test segment of a wind tunnel. In one experiment, the total frictional force exerted by the wind on the plate was measured,and the shear stress (*) was obtained by calculating the total frictional force per unit of plate area,and simultaneously measuring the wind velocity (U) at an elevation of =52cm from the bottom of the wind tunnel, the air pressure, and the air temperature.The air density was calculated using the Clapeyron equation, and the shear velocity was calculated using u*=(*/g)1/2. In these experiments, the densities of the roughness elements were 0, 25, 100, 200, 400, 800, 1300, 1600, and 2400 m-2, and the wind velocities (U) were5, 6, 7, 8, 9, 10, 11, 12, and 13ms-1. Dong et al. ignored the displacement height(d) in equation (3), and calculated the roughness length (z0) from experimental data for u* andU. The results are shown in Figure 2 with the exception of data for n=0 and n=25.
Fig.2 shows that z0 decreases with increasing wind speed, butu*proportionally increases with increasing wind speed, and that each of z0 and u*increases at a progressively slower rate as n increases. We did not find that z0 or u* decreased with increasing n, so it appears that the isolatedroughness flow and the wake interference flow, but not the skimming flow, were included in the results reported by Dong et al. However, we did observe an important dynamical behavior: that the increase in z0 and u*with increasing n became progressively slower and slower as n increased.
(a) (b)
Fig.2. Roughness length (a) and shear velocity (b) as a function of the wind velocity (U) above the turbulent boundary layer based on empirical data.
Equation (20) can be rewritten as follows
（20）
We fittedthe experimental data for each wind velocity to equation (20), and list the results in Table 1.
Table 1 Cr0 and rvalues in equation (20)
U(m/s) / 5 / 6 / 7 / 8 / 9 / 10 / 11 / 12 / 13Cr0
r / 0.9594
2.1296 / 0.8940
2.1217 / 0.8473
2.1187 / 0.7911
2.1070 / 0.7486
2.0822 / 0.7113
2.0730 / 0.6693
2.0320 / 0.6432
2.0043 / 0.6220
2.0013
In Table 1, we see that Cr0 decreases continuously with increasing U, (i.e., Cr0 decreases withincreasing Re=Uh/v ), but that r remains nearly constant (i.e., only decreases by around 6%) for the range of wind velocities that were studied, and this suggests that equation (20) successfully describes the dynamics of roughness and shear velocity.However, one clearunsuccessful aspect of equation (20) is shown in Fig.2:the expected decrease in roughness length and shear velocity doesn’t appear until n=2400, butr-1=0.4821 is obtained by fitting the experimental data to equation (20), and this suggests that roughness length and shear velocity should decrease with increasing distribution density of the roughness elements when n1004. Here, we must emphasize an important experimental characteristic, namely that the rate of increase of z0 and u* with increasing n becomes too slow when n800.
Based on the above analysis, we propose a hydro-mechanical model for the roughness arrays created from thin cylindrical sticks with different distribution densities (n) arranged on a smooth wall.
6. Dynamical nature of roughness length and shear velocity with thin cylindrical roughness elements
6.1 Momentum equations for the different flow sub-layers
Suppose that thin cylindrical elements, with Eh(where E is the diameter of the cylindrical elements and his their height), are arranged on a smooth wall, and that the volume fraction of the roughness elements in one unit volume is near zero. This situation is comparable to many real cases; for example, the upper limit of the volume fraction was only 1.7% in Dong et al.'s experiment[5]for different distribution densities of the roughness elements.
We can define a new shear velocity function u*(z)
（21）
（22）
Thus, we have
（23）
Where=1-nπE2/4. Similarly, we can define a new Reynolds number function
（24）
The results of many experiments have shown the existence of three regions based on the value of the Reynolds number: Re*(z) ≤ 5, 5Re*(z) ≤ 30 and Re*(z)30, respectively, represent the laminar, transitional, and turbulent sub-layers. Let us ignore the transitional sub-layer, and assume that the laminar sub-layer changes directly into the turbulent sub-layer once Re*(z)5. Defining the root of equation Re*(z)=5 as a characteristic height,h*, and then, using Prandtl's mixing-lengththeory, we can divide the boundary layer into different sub-layers as follows.
When h* ≥ , the boundary layer can be described as follows
（25）
When hh*, the boundary layer can be described as follows
（26）
When h*=h, the boundary layer can be described as follows
（27）
And when h*h, the boundary layer can be described as follows
（28）