In PHOENICS the effect of roughness is taken into account by allowing the roughness parameter E, in the logarithmic law for the near-wall velocity, to vary according to the empirical laws proposed by Jayatilleke [1969] for sand-grain roughness.
This enhancement is available for use with both equilibrium and non- equilibrium logarithmic wall functions, whether applied at wall-type patches, or at solid-fluid interfaces if EGWF is set to T in Q1 (see EARTH_Generated_Wall_Functions entries in PHENC).
For a smooth wall, E is a constant, denoted Em, taken as 8.6; and for a rough wall E is expressed as a function of the roughness Reynolds number, Rer, defined as:
Rer = Uτ*hr/νl
(8.4.1)
in which hr is the absolute value of the equivalent sand-grain roughness height.
The formula for E is as follows:
when Rer <3.7: E=Em;
when 3.7 <Rer < 100: E=1./√(a*(Rer/b)2+(1-a)/Em2);
and when Rer > 100: E = b/Rer
where: b = 29.7; a = (1.+2*X3-3*X2); and X = 0.02248(100.-Rer)/Rer0.564
(8.4.2)
Allowance is made for rough walls in the calculation of the Stanton number by replacing the smooth-wall sublayer resistance function, Pm, with the empirical formula of Jayatilleke [1969]:
P = 3.15σl0.695(1/E-1/Em)0.359+Pm(E/Em)0.6
(8.4.3)
where Pm is the value of P under fully smooth conditions.
The sand-grain roughness height is the roughness height which gives the same actual resistance coefficient as that caused by the real non-uniform wall roughness. The height of this equivalent sand-grain roughness is often specified in terms of a relative roughness (e.g. equivalent protrusion height-to-diameter ratio in pipe flows).
Further information giving recommendations and formulae for the effective sand-grain roughness height can be found in the literature (see for example H.Schlichting(1968)).
PHOENICS provides a set of fully-rough wall functions that are suitable for a near-wall layer in local equilibrium defined in terms of the effective aerodynamic roughness height y0 and the zero-plane displacement height d. The most common application area for this type of wall function is the atmospheric boundary layer under neutral conditions.
These wall functions are similar to the "sand-grain-roughness" equilibrium wall functions described earlier, except that they now take the form:
Ur=(Uτ/κ)*ln[(y-d)/y0]
k=Uτ2/√ℂμ
ε=Uτ3/[κ*(y-d)]
(8.4.4)
(8.4.5)
(8.4.6)
where Ur is the resultant velocity parallel to the surface, Uτ is the friction velocity, κ=0.41 is von Karman's constant; ℂμ=CμCd=0.09; y0 is the effective roughness height, and d is the zero-plane displacement height due to vegetation, buildings or other obstacles.
The roughness height y0 is related to the size of the roughness elements on the surface, and is typically between 1/10 and 1/30 of the average height of the roughness elements. Some commonly-used values of y0 are listed below:
| Surface type | Roughness height yo (m) |
| Calm open water | 0.0002 | Rough open sea | 0.001 |
| Open flat terrain, grass, few isolated obstacles | 0.03 |
| Low crops, occasional large obstacles | 0.10 |
| High crops, scattered obstacles | 0.25 |
| Parkland, bushes, numerous obstacles | 0.50 |
| Suburb, forest, regular large obstacle coverage | 0.50 to 1.0 |
Values greater than 1m are rare and indicate excessively rough terrain.
The zero-plane displacement d is the height above the ground at which zero wind speed is achieved as a result of flow obstacles such as trees or buildings, i.e. Ur=0.0 at y=y0+d. Often the displacement height is zero, but for flow over an array of densely packed objects (e.g. a forest, cropland or buildings), an offset in height is introduced into the log law to allow for the upward displacement of the flow by the surface objects. This displacement height d is usually estimated as 2/3 of the average height of the obstacles.
The displacement height d is a positive quantity, although the user can set d=-y0 to use wall functions that are consistent with the inlet wind profiles proposed by Richards & Hoxey (1993), which are often used in wind-engineering simulations.
Equation (8.4.4) implies that the absolute value of the wall shear stress can be computed from:
τw = sρUr2
where
s = (κ/ln[(y-d)/y0]2
(8.4.7)
(8.4.8)
If x is a coordinate direction aligned with the surface, the stress in that direction is:
τw,x = sρUrUx
(8.4.9)
where Ux is the x-direction velocity at the near-wall grid node at distance y from the wall.
For heat and mass transfer, the Stanton number is computed from a modified Reynolds analogy, i.e. St=s/σt.
In PHOENICS there are two methods by which the user can set the roughness height for
walls:
This roughness height is then applied to all walls in the solution domain (whether defined by wall-type patches, or found by EARTH if e.g.WF=T), unless the roughness for a particular wall is defined by the alternative method described below.
The default value of WALLA is 0.0, which defines a hydrodynamically smooth wall.
PATCH(NAME,WWALL,.,.,.,.,.,.,.,.)
SPEDAT(SET,ROUGHNESS,NAME,R,0.001)
COVAL(NAME,U1,GRND2,0.0); etc.
This allows the user to modify the roughness height in GROUND coding; and, for example, introduce an arbitrary distribution of roughness height across a patch. The appropriate subroutine can be called from any convenient place in GROUND, except Group 1 Section 1.
For details of the work with EARTH patch-wise variables see entry PATCH-WISE VARIABLES entry to PHENC.
The displacement height is set through the PIL variable WALLB stored in the common block /RSG/ of the GRDLOC include file.
Note that the coding for the wall functions can be found in the file GXWALL which comprises the subroutines GXWFUN, FNSKIN, FNCOEF and GXYPLS.