1 LOW-MOBILITY TRANSPORT OF COARSE-GRAINED BED MATERIAL UNDER WAVES AND CURRENTS

Rubble mound near-bed structures are among the more common hydraulic structures. Examples include bottom protections to prevent scour near bridge piers, offshore structures, weirs and sluices and at the toe of bank protections or breakwaters, but also pipe covers on offshore pipelines. The cornerstone of the theory behind the design of these structures was laid by Izbash (1930) and, most notably, Shields (1936), who linked the stability of a stone in the near-bed structure to the shear stress exerted on the bottom by a flow. Even today, some 75 years later, the design of near-bed structures is still largely based on Shields' work. Shields postulated that for values of his stability parameter below a certain critical value (the socalled 'threshold of motion') the stones in the bed would not move at all. In reality stones do not exhibit this kind of behaviour; stones have been shown to move at any value of the Shields parameter, also below the 'threshold of motion'. The mobility of the stones does increase with higher values of the Shields parameter, indicating that this is really a mobility parameter rather than a stability parameter. In practice, this concept has given rise to a conservative ‘statically stable’ design approach: a threshold of motion is selected and the size and weight of the stones are chosen such that the Shields parameter does not exceed the critical value. Alternative, ‘dynamically stable’ design approaches in which some movement of stones is allowed in combination with an appropriate maintenance program may be more economical. However, these are not commonly adopted because the transport of the bed protection material cannot be accurately predicted and thus the rate of deterioration of the bed protection and the required intensity of the maintenance programme cannot be assessed. Research has been undertaken to investigate the transport of bed protection material and to improve the feasibility of a ‘dynamically stable’ design approach. Application in marine environments requires that the loading consists of waves and currents, or a combination of these.


INTRODUCTION
Rubble mound near-bed structures are among the more common hydraulic structures.Examples include bottom protections to prevent scour near bridge piers, offshore structures, weirs and sluices and at the toe of bank protections or breakwaters, but also pipe covers on offshore pipelines.
The cornerstone of the theory behind the design of these structures was laid by Izbash (1930) and, most notably, Shields (1936), who linked the stability of a stone in the near-bed structure to the shear stress exerted on the bottom by a flow.Even today, some 75 years later, the design of near-bed structures is still largely based on Shields' work.
Shields postulated that for values of his stability parameter below a certain critical value (the socalled 'threshold of motion') the stones in the bed would not move at all.In reality stones do not exhibit this kind of behaviour; stones have been shown to move at any value of the Shields parameter, also below the 'threshold of motion'.The mobility of the stones does increase with higher values of the Shields parameter, indicating that this is really a mobility parameter rather than a stability parameter.
In practice, this concept has given rise to a conservative 'statically stable' design approach: a threshold of motion is selected and the size and weight of the stones are chosen such that the Shields parameter does not exceed the critical value.Alternative, 'dynamically stable' design approaches in which some movement of stones is allowed in combination with an appropriate maintenance program may be more economical.However, these are not commonly adopted because the transport of the bed protection material cannot be accurately predicted and thus the rate of deterioration of the bed protection and the required intensity of the maintenance programme cannot be assessed.
Research has been undertaken to investigate the transport of bed protection material and to improve the feasibility of a 'dynamically stable' design approach.Application in marine environments requires that the loading consists of waves and currents, or a combination of these.

PREVIOUS WORK Paintal formula
Transport of rock is characterized by low mobility: there is occasional movement of a stone, rolling or sliding over a short distance.Also, for bed protections the individual movement of a single stone, or a few stones, can determine the functional stability of a structure.This requires a different level of analysis compared to most available morphological transport prediction formulas that focus on transport of sand.In addition, the difference in physical phenomena like bed roughness, grain shape (roundness), turbulence and porous flow through the bed make available sand transport formulas inapplicable for analysis of rock transport.Paintal (1969) conducted measurements on just these issues (coarse-grained, low-mobility transport) for steady currents; his dataset is still one of the few available in this field.Paintal developed a theoretical probabilistic model of stone transport in which stones move in discrete 'steps', rolling and sliding a short distance along the bottom before being redeposited again.The distribution of these steps is random in time and space, and there is no general ongoing movement.Paintal found from this model that the dimensionless transport rate is a function of the Shields parameter only (Paintal 1969).From laboratory measurements Paintal then fitted the transport equation: where  q = q s /√(gΔd n50 ) is the dimensionless transport rate, Ψ = u * 2 /(gΔd n50 ) is the dimensionless bed shear stress (Shields parameter), q s is the transport rate through a cross section per unit width (in m 3 /m/s), Δ = (ρ stone -ρ water )/ρ water is the relative density of the stones, d n50 is the nominal median diameter of the stones (in m) and u * is the bed shear stress velocity (in m/s).
This formula was later corrected by WL|Delft Hydraulics (Mosselman and Akkerman 1998) for temperature effects and the apparently false assumption of a hydraulically rough bottom by Paintal in one of his test series.The corrected formula yields: So the second formula remains unaltered, only the transition point ( = 0.085) has moved to ensure a correct transition between the two formulas.Perhaps confusingly, the Paintal formula was corrected a second time by Mosselman and Akkerman (1998).This second correction yields the following formula: In this research we have further used this last version of the Paintal formula.
As said, these formulas have not yet been very commonly applied in engineering practice.Looking at the exponent of  in these formulas (and especially the original Paintal formula) may explain why: the transport rate varies with the 16th power of the shear stress, and since   u 2 , it varies with the 32nd power of the velocity.Apart from the strange, perhaps unfamiliar and 'suspicious' magnitude of this exponent (more classical morphological formulas, for higher -values, predict transport rates varying with the 5th power or so to the velocity) this also means that a small inaccuracy in the estimate of the velocity results in large deviations of the predicted transport.Paintal himself seems to be aware of this consequence of his formula when he writes that "the 16th power correlation […] stipulates that a small increase in shear stress causes a noticeable change in bed load transport" (Paintal, 1969)

Additional validation measurements
These matters have been investigated further by Royal Boskalis Westminster and Delft University of Technology (DUT) in the late 1990s.De Boer (1998) conducted additional transport measurements for low hydraulic loads (still limiting the application to steady currents) in the DUT laboratory facilities and found much lower transport values than Paintal did.The differences between Paintal and De Boer have, at least partly, been explained by Forschelen (1999), who found that they could to a large extent be attributed to a different definition of bottom roughness by the two researchers and to an underestimation of friction effects caused by the side walls of the flume in De Boer's experiments.Forschelen also did some transport measurements of his own, which seem to be more in line with Paintal's original findings.The data of Paintal, De Boer, Forschelen and equations (1) through (3) above are illustrated in figure 1.

Wave-current interaction
The above methods are all based on the bed shear stress as the load parameter.When attempting to extend the applicability to a combination of waves and a current, the combined wave-current bed shear stress must be determined.This requires consideration of the following: 1. Spatial orientation: the waves and current may propagate at an angle, which implies a vector analysis of the combined shear stress.One could analyze the combined shear stress "in the direction of the flow", taking the projection of the combined shear stress vector in the flow direction, or one could use the absolute values of the combined shear stress.Given the limitations of the facilities used in this present research we limit ourselves to collinear flow (φ = 0 o during half the wave cycle, and φ = 180 o during the other half).The difference between analysis "in the direction of the flow" and absolute values is then still relevant in case the combined shear stress is directed against the flow direction during some part of the wave cycle.2. The cyclic nature of the wave shear stress: one could analyze either the maximum (amplitude) of the combined shear stress or an average value over the wave cycle (or alternatively over half the wave cycle) 3. Wave-current interaction (WCI): Soulsby (1997) lists three ways in which the waves and the current interact: (a) modification of the phase speed and wavelength of the waves by the current, leading to refraction of the waves, (b) interaction of the wave and current boundary layers, leading to enhancement of both the steady and oscillatory components of the bed shear stress and (c) generation of currents by the waves, including longshore currents, undertow and mass transport (streaming) currents.In this present paper we focus on the second phenomenon only.
Soulsby et al compared these eight models in terms of their prediction of the mean and maximum combined shear stress for the same input parameters.It can be deduced that 'mean' in the context of item (1) and (2) above should be interpreted as "average over the whole wave cycle in the direction of the current" (Van den Bos 2006).The input parameters are the same for all models: the relative bottom roughness z 0 /h, the relative wave excursion a 0 /z 0 and the angle between the current and the wave propagation direction .Soulsby et al give their results in plots of the dimensionless input parameters X = τ c / (τ c + τ w,max ) against the dimensionless mean combined shear stress Y = τ wc, mean / (τ c + τ w,max ) and dimensionless maximum shear stress Z = τ wc, max / (τ c + τ w,max ).In this notation τ c is the current-only bed shear stress, τ w,max is the amplitude of the wave-only bed shear stress and τ wc is the combined wave-current bed shear stress.The parameter X can be seen as a measure for the relative strength of the current and wave components, for X = 0 there are only waves, for X = 1 there is only a current.A typical plot from the comparison in Soulsby et al is given in Figure 2. As can be seen from this plot, there are quite some differences in the predictions from the various models.It appears that, roughly speaking, three groups of models can be discerned (Bijman, 2000) 1. the first group, containing Bijker, and Van Kesteren and Bakker, give strikingly larger maximum shear stresses than the other models, and also higher mean stresses than the other models for wave dominated situations (low X values); 2. The second group is formed by Grant and Madsen, Myrhaug and Slaattelid, and Christofferson and Jonsson, who predict higher values than the other models for current dominated situations; 3. Finally, Fredsøe, Davies et al and Huynh-Thanh and Temperville predict lower values than the other models for all situations (and their predictions appear to match each other closely).This present research has investigated which WCI model is best suited to predict stone transport.For simplicity reasons only a few WCI models have been applied in the comparison, being those by Bijker (1967), Grant and Madsen (1979) and Fredsøe (1984).The choice is mainly based on the above distinction into three groups and one representative model from each group has been selected.Also a more physical reasoning has been applied to support this selection, as described below.
The existence of the wave boundary layer, and the extra turbulence it generates, is felt by the current as increased resistance which explains the enhanced current-only bed shear stress.On the same time the existence of a steady current influences the growth of the wave boundary layer.The dominating length scales and the turbulence characteristics are different for uniform currents and for waves; this is one of the main reasons why simply adding the resulting shear stresses is not a valid approach.A physically more correct approach would be to assume a combined wave-current wave boundary layer and treat the flow inside this layer as a separate flow problem, with separate boundary conditions and turbulence closure.This is exactly what the various wave-current interaction models seek to do.The differences between the models are mainly caused by the way in which they incorporate turbulence inside the boundary layer.This provides us with a way to categorize the wavecurrent interaction models (Van Rijn 1993).
All wave-current interaction models discussed here use the eddy viscosity concept in their turbulence closure.Generally there are four ways to apply this concept which are, in order of increasing complexity: 1. assuming a constant eddy viscosity, which is equivalent to using the mixing length hypothesis and assuming that l m is constant.This leads to a linear velocity profile in the wave boundary layer.This is the basis of the Bijker model.2. using the mixing length hypothesis and assuming that l m varies with the distance from the bottom, which leads to a logarithmic velocity profile in the boundary layer.Both Grant and Madsen, and Fredsøe use this concept in their models, the main difference being that Grant and Madsen assume the boundary layer thickness and eddy viscosity to be constant in time, while Fredsøe allows these parameters to vary with the wave cycle.3. using a balance equation for the turbulent kinetic energy and an estimated length scale (oneequation model).This approach was adopted by Davies, Soulsby and King.4. using a balance equation for the turbulent kinetic energy and an (assumed) balance equation for the dissipation rate of the turbulent kinetic energy (k- closure), as was done by Huynh-Thanh and Temperville.
So we see that each of the selected WCI models represents another turbulence closure method, with increasing complexity.The models using numerical turbulence closure methods (items 3 and 4 in the above list) are omitted from the present research as their predictions do differ much from the Fredsøe method.This increasing complexity is also found in the practical application of the models: the Bijker model results in an explicit straightforward equation, the Grant and Madsen model requires iteration and the Fredsøe model requires the (numerical) solution of a rather complicated differential equation within an iteration loop; a more elaborate overview is given in Van den Bos (2006).
MODEL TESTS Bijman (2000) conducted transport measurements using colored strips of bed material in the facilities of Delft University of Technology (length 45 m, height 1.0 m, width 0.80 m).The experiments were carried out using steady current and regular waves, and water-working has been applied (1800 seconds uniform flow followed by 2400 waves prior to testing) to avoid bias of the test results by accidentally poorly positioned stones on the bed.The test setup is described in Figure 3. Bijman uses two measurement zones, each containing ten 5 cm wide strips of loose colored stones with a layer thickness of 5 cm.A strip of 15 cm along both flume walls was not taken into account in the analysis, so the effective width of the measurement areas was 50 cm.The material in between the measurement zones is glued to the bed.Transport was measured by counting the individual stones that were moved outside their original strip and expressed as the volume per unit width and time q s (i.e.divided by the effective width of the measurement zone and the test duration) and eventually as the dimensionless transport parameter  q (as introduced above).
The range of parameters that were varied in his research is given in table 1.The wave Reynolds number is defined as Re w = û 0 a 0 / in which û 0 is the amplitude of the horizontal bed orbital velocity (m/s), a 0 is the amplitude of the horizontal particle excursion at the bed (m) and  is the kinematic viscosity of the fluid (m 2 /s).This parameter can be used in combination with the relative bed roughness a 0 /k s to classify the flow regime of these tests (Jonsson 1966), see figure 4.
In nature, the flow can be assumed to be rough turbulent; it turns out that Bijman's experiments were carried out at the very edge of the rough turbulent zone and in the transition zone from laminar to rough turbulent.Therefore viscous scale effects may play a role in his experiments to some extent.

Remarks
The bed roughness is taken as The combined shear stress in this overview was calculated with the Bijker method: The particle Reynolds number is based on the combined shear stress velocity The current dominance parameter X is defined according to Soulsby

ANALYSIS
As a first step in the analysis we have determined the general predictive power of the available methods without looking at the Paintal formula itself.That is, we have fitted the Bijman data to a generalized transport equation  q = a . wc b in which the combined wave-current dimensionless shear stress  wc is calculated using the three selected WCI models and four combinations of the possible spatial and temporal averaging techniques described above.In addition we have made calculations without wave-current interaction as a baseline for comparison.In this first stage we are not interested in the values of the regression coefficients a and b; instead the focus is on the statistical goodness-of-fit r 2 (coefficient of determination).This parameter is used as a measure for the general ability of a WCI model / averaging technique combination to predict the observed transport rates.The analysis will then be continued with the highest-scoring combinations.Bryman and Cramer (2005) give as a guideline that the goodness-of-fit can be classified as 'modest' for r 2 = 0.16 -0.49(r = 0.40 -0.69), 'high' for r 2 = 0.49 -0.81 (r = 0.70 -0.89) and 'very high' for r 2 > 0.81 (r > 0.90).The results are given in table 2. This table shows that 1.The goodness-of-fit increases considerably when wave-current interaction is taken into account, which confirms the theoretical finding that this phenomenon cannot be ignored.2. The more advanced wave-current interaction models by Grant and Madsen, and Fredsøe, fit the data better than the simpler model by Bijker, although the difference (especially in the combinations with the average shear stress in the current direction) is not large.3. The Fredsøe model performs just as well as the model by Grant and Madsen when the average shear stress in the current direction is used; it performs better for all other averaging procedures.4. For all shear stress models, the stability parameters based on the maximum shear stress perform worse than stability parameters based on the average shear stress.5.For the Fredsøe model the average shear stress in the current direction performs just as good as the other two averaging procedures, for all other models the average shear stress in the current direction provides the best fit.
From this first step it is concluded that the average shear stress in the direction of the current is the best predictor of stone transport.This is the same parameter as the 'mean' shear stress as analyzed by Soulsby (1997).In terms of WCI models it is concluded that the prediction improved when WCI is taken into account, and all the three analyzed models roughly perform similarly.The analysis is continued with these three models and only the average shear stress in the current direction.The other averaging techniques are discarded.In the second step the predictions from the remaining models are compared to the Paintal formula (in the Mosselman and Akkerman correction, equation 3).The reason is that from a practical engineering point of view it would be desirable to have one single stone transport prediction formula for both current-only situations and wave-current combinations.Since the applicability of the Paintal formula for the first case is not disputed in this paper, it is preferable to find a wave-current prediction model that results in the same formula over having a second formula with different coefficients, even if this second formula has a (slightly) higher statistical goodness-of-fit.
To this purpose, the Bijman data were plotted in a ( wc ,  q ) graph along with the Paintal formula, with  wc calculated for each of the three remaining methods, see figure 5.

CONCLUSION
From the above analysis it is concluded that the Paintal formula (in the version by Mosselman and Akkerman (1998)) that was developed to predict the low-mobility transport of coarse-grained bed material under steady current load, can also be used to predict the transport under a combination of waves and a steady current, if the resulting bed combined shear stress is calculated as the mean shear stress (average shear stress in the direction of the current) using the wave-current interaction model of Fredsøe (1984).This conclusion is valid within the limitations of the underlying dataset, see table 1.

DISCUSSION
The applicability of the Paintal formula to predict transport of stones from a bed protection has been extended to loading by waves and currents.The data cloud analyzed in this paper is presented in the context of the previous work and the prediction formulas in figure 6.It was concluded that the best predictor for stone transport was the average combined shear stress in the direction of the current.From a physical point of view this predictor can be seen as the enhanced current shear stress, i.e. the shear stress exerted by the current but enhanced by the presence of the waves (by means of their interference with the current boundary layer).This gives a physical explanation to the observed phenomena: the current is still responsible for the transport of the stones, the presence of the waves only acts as a catalyst.
The best predictor of wave-current interaction in combination with the Paintal formula was found to be the method proposed by Fredsøe (1984).This is the most complex of the methods analyzed in this research.In engineering practice this method has limited applicability since it requires the numerical solution of a differential equation within an iteration loop.Therefore for practical applications the parameterization method developed by Soulsby et al (1993) is recommended.This method gives an explicit formula that approximates the results from the full model.
Finally, it is noted that the methods discussed in this paper are all based on the bed shear stress as the acting load on the stones.Other methods, using for instance velocities, accelerations or direct turbulence parameters also exist.In this present paper we follow the engineering practice in which the shear stress is commonly used.A discussion of the relative merits and disadvantages of alternative methods is considered outside the scope of this present paper.A more extensive overview of these methods is given in the full report on which this paper is based (Van den Bos 2006).

Figure 1 .
Figure 1.Illustration of datasets for low mobility transport and fitted equations, for steady currents.Note that original Paintal data points (open squares) are uncorrected values which do not fit the (corrected) Mosselman and Akkerman formula

Figure 2 .
Figure 2. Example of prediction of maximum and mean combined shear stress for various wave-current interaction models (Soulsby et al 1993)

Figure 4 .
Figure 4. Model test flow regime

FredsoeFigure 5 .
Figure 5. Bijman (2000) data set plotted against the Paintal formula for three WCI models using the average combined shear stress in the current direction q

Figure 6 .
Figure 6.Extension of the applicability of the Paintal formula to loading by waves and currents