Improvement of the Full-Range Equation for Wave Boundary Layer Thickness

In order to improve the accuracy of the original full-range equation for wave boundary layer thickness, with special reference to increasing its applicability to tsunami-scale waves, a theoretical investigation is carried out to derive a dimensionless expression which is valid under both smooth and rough turbulent regimes. A coefficient in the equation is determined through a comparison with k-ω  model computation results for tsunami-waves along with laboratory scale oscillatory flow experiments. Thus, the improved full-range equation for wave boundary layer thickness enables us to cover a wide range of wave periods from wind-wave to tsunami.


Introduction
The 2011 Tohoku tsunami caused highly devastating disasters to human lives, infrastructure, and coastlines in the north-eastern region of Japan. The extreme tsunami-induced flow velocity caused remarkable morphology change in coastal and estuarine areas in Tohoku region. Using highly valuable field data of morphology changes [1,2], many researches have been able to calibrate and validate a numerical model for tsunami-induced sediment transport and resulting morphology change [3,4]. In these numerical investigations, it is common to use a steady friction factor for estimating bottom shear stress under tsunami. However, according to a numerical study using k-ω turbulence model by Williams and Fuhrman [5], tsunami-induced boundary layer is wave-like in the sense that the boundary layer does not necessarily span the entirety of the water depth. In addition, a recent field investigation during the 2010 Chilean tsunami revealed that the velocity profile at a water depth of 10 m is extremely similar to that under wave motion with very thin boundary layer thickness [6]. Tinh and Tanaka [7] reported that the ratio of δ/h (δ: the wave boundary layer thickness under tsunami, h: the water depth) is suitable for judging whether a tsunami-induced boundary layer is wave-like or current-like.
Wave boundary layer thickness has been formulated mainly using experimental results obtained in an oscillating tunnel [8][9][10][11]. It is well known that laminar wave boundary layer thickness is proportional to (νT) 1/2 (ν: kinematic viscosity, T: wave period) according to Stokes second problem [12]. Thus, the boundary layer thickness developed in a wave flume is rather small (order of 0.1 cm). For this reason, a closed type of oscillating tunnel has been developed [13] to increase the period of oscillatory motion, and has been extensively applied by many researchers for studying various topics related to wave motion including boundary layer structure, sediment movement, and force acting on coastal structures under oscillatory motion [14][15][16][17].
Nonetheless, the period of the oscillatory motion in these experimental facilities is extremely smaller than that under tsunami. Tinh and Tanaka [7] carried out an investigation on the development of a bottom boundary layer under shoaling tsunami using the formula proposed by Sana and Tanaka [11]. In their study, a hypothetical tsunami shoaling case was configured. Similar to Williams and Fuhrman [5], an initial offshore tsunami of 1 m high with the wave period 15 min at the water depth of 4000 m that propagates over 200 km to the shallow area is considered. As a result, it is shown that depth-limited condition, δ = h is not satisfied by tsunami boundary layer even at the water depth shallower than h = 10 m, which shows good agreement with the field measurement by Lacy and colleagues [6]. However, it is noted that the equation proposed by Sana and Tanaka [11] is applicable to wind-generated waves, and it is not clear whether it can be applied to long period wave by extrapolation. In the present study, therefore, a theoretical framework is developed for boundary thickness under turbulent conditions covering both wind-waves and tsunami. Furthermore, a semi-empirical formula is proposed for wave boundary layer thickness using laboratory data along with numerical results by k-ω turbulence model.

Comparison of Existing Formulae
Our interest herein will be limited to the wave boundary layer under sinusoidal wave motion. In the case of the laminar boundary layer, the governing equation is given as follows.
where u: the horizontal velocity in the wave boundary layer, t: the time, U: the free stream velocity outside the boundary layer, and z: the vertical coordinate. In order to solve Equation (1), the following two boundary conditions are imposed.
where U m : the amplitude of U, and σ: the angular frequency (=2π/T). From these equations, the following analytical solution can be obtained for the velocity distribution in a boundary layer [16].
where δ 1 : the Stokes layer thickness (= √ 2ν/σ). The time-variation of velocity distribution over a wave cycle from Equation (4) is shown in Figure 1. One of the important characteristics of the wave boundary layer is an overshoot phenomenon near the bottom as illustrated in Figure 1. Jensen et al. [9] defined the height of overshoot under the wave crest as a wave boundary layer thickness δ, and obtained the following analytical expression.
where the subscript (L) denotes laminar condition, and Re in Equation (5) is Reynolds number defined by in which a m is the amplitude of water particle movement outside the wave boundary layer as At higher Reynolds numbers, transition to a turbulence occurs around Re = 2.5 × 10 5 [19,20]. Depending on different criteria for defining the transition, slightly different values have been proposed (see [21,22]). For turbulent regime, analytical solution cannot be obtained for boundary layer thickness. Hence, an empirical approach has been applied based on laboratory experiments or numerical analysis using a turbulence model such as k-model and k-model. The following equations have been proposed for smooth turbulent and rough turbulent conditions. Alternative definitions for wave boundary layer thickness have been employed by other researchers, such as the minimum distance from the bottom to a point where the horizontal velocity equals free-stream wave velocity amplitude [16], and the elevation of 1% or 5% velocity deficit at the outer edge of the boundary layer [8,18]. Among these definitions, that of Jensen et al. [9] is the most suitable for the present study, since we will discuss the depth-limited condition of tsunami boundary layer, in which a shear-free condition needs to be exactly satisfied on the water surface.
According to Equation (5), it is easily shown that wave boundary layer thickness under laminar condition is proportional to (νT) 1/2 as mentioned before.
At higher Reynolds numbers, transition to a turbulence occurs around Re = 2.5 × 10 5 [19,20]. Depending on different criteria for defining the transition, slightly different values have been proposed (see [21,22]). For turbulent regime, analytical solution cannot be obtained for boundary layer thickness. Hence, an empirical approach has been applied based on laboratory experiments or numerical analysis using a turbulence model such as k-ε model and k-ω model. The following equations have been proposed for smooth turbulent and rough turbulent conditions.  (13) where the subscripts (S) and (R) denote rough turbulent and smooth turbulent, respectively, and k s is the equivalent roughness. Among these, Williams and Fuhrman's [5] equations have been proposed based on k-ω numerical simulation for tsunami-like long period waves. Figure 2 shows comparison among Equations (8)-(10) for smooth turbulent condition. It is seen that Equations (8) and (10) have similar behavior, whereas the equation by Sana and Tanaka [11] has a smaller dependence on the Reynolds number. The line designated "Present study" represents a new equation proposed in this paper, which will be described later.  (13) where the subscripts (S) and (R) denote rough turbulent and smooth turbulent, respectively, and ks is the equivalent roughness. Among these, Williams and Fuhrman's [5] equations have been proposed based on k-numerical simulation for tsunami-like long period waves. Figure 2 shows comparison among Equations (8)-(10) for smooth turbulent condition. It is seen that Equations (8) and (10) have similar behavior, whereas the equation by Sana and Tanaka [11] has a smaller dependence on the Reynolds number. The line designated "Present study" represents a new equation proposed in this paper, which will be described later. Similarly, Figure 3 shows comparison of three equations for boundary layer thickness under rough turbulent condition. Again, the equation proposed by Sana and Tanaka [11] differs distinctly from other equations especially at higher / due to its larger power in Equation (12), which may cause underestimation of the boundary layer thickness for tsunami waves with a longer wave period. Similarly, Figure 3 shows comparison of three equations for boundary layer thickness under rough turbulent condition. Again, the equation proposed by Sana and Tanaka [11] differs distinctly from other equations especially at higher a m /k s due to its larger power in Equation (12), which may cause underestimation of the boundary layer thickness for tsunami waves with a longer wave period.
It must be stressed here that each empirical equation is applicable in the region of available data on which the equation is based. For this reason, applicability of equations needs to be carefully investigated when it is extrapolated to longer wave motion like tsunami. It must be stressed here that each empirical equation is applicable in the region of available data on which the equation is based. For this reason, applicability of equations needs to be carefully investigated when it is extrapolated to longer wave motion like tsunami.

Theoretical Consideration
Hereafter, dimensionless expression will be derived for boundary layer thickness which is valid for both smooth turbulent and rough turbulent conditions.
As described above, the boundary layer thickness under laminar condition can be expressed as ∝ ( ) 1/2 . This is a representative length common for various physical processes expressed in terms of the diffusion equation. For example, equivalent representative length can be derived in heat conduction theories [23], viscous fluid theories [12] and also in a shoreline change model based on a one-line model [24,25]. After transition to turbulence at higher Reynolds numbers, the molecular viscosity in the equation for wave boundary layer thickness might be replaced by eddy viscosity . Kajiura [26] and Tanaka and Shuto [27] used the following assumption to derive an analytical solution for turbulent wave boundary layer.
where : the Karman constant (=0.4), and * : the maximum friction velocity. The dependence on in Equation (14) is originated from mixing length, = , due to turbulent momentum transfer. Alternative assumption for wave-induced eddy viscosity may be constant expression, rather than spatially varying as Equation (14), as assumed by Kajiura [26] in the outer layer of the wave boundary layer. Then, replacing in Equation (14) by , the following equation can easily be obtained for wave boundary layer thickness under turbulent conditions. * = .
where * is defined by

Theoretical Consideration
Hereafter, dimensionless expression will be derived for boundary layer thickness which is valid for both smooth turbulent and rough turbulent conditions.
As described above, the boundary layer thickness under laminar condition can be expressed as δ ∝ (νT) 1/2 . This is a representative length common for various physical processes expressed in terms of the diffusion equation. For example, equivalent representative length can be derived in heat conduction theories [23], viscous fluid theories [12] and also in a shoreline change model based on a one-line model [24,25]. After transition to turbulence at higher Reynolds numbers, the molecular viscosity in the equation for wave boundary layer thickness might be replaced by eddy viscosity ν t . Kajiura [26] and Tanaka and Shuto [27] used the following assumption to derive an analytical solution for turbulent wave boundary layer.
where κ: the Karman constant (=0.4), and u * m : the maximum friction velocity. The dependence on z in Equation (14) is originated from mixing length, l = κz, due to turbulent momentum transfer. Alternative assumption for wave-induced eddy viscosity may be constant expression, rather than spatially varying as Equation (14), as assumed by Kajiura [26] in the outer layer of the wave boundary layer. Then, replacing z in Equation (14) by δ, the following equation can easily be obtained for wave boundary layer thickness under turbulent conditions.
where a * m is defined by where f w : the wave friction coefficient. Here, the relationship between a friction coefficient and friction velocity u * m is given by the following formula.
It is interesting to note that the non-dimensional quantity defined by Equation (15) has already been used in theoretical study for turbulent wave boundary layer by Kajiura [26], Grant and Madsen [28], and Tanaka [29]. Another interesting point is mathematical similarity between laminar and turbulent flow in the denominator of the left-hand side of Equation (5), a m , and the corresponding expression in Equation (15), a * m .
In order to validate Equation (15), velocity profile data under wave crest (σt = 0 deg) were collected as summarized in Table 1. In both laboratory experiments by Sleath [8] and Jensen [30], a U-shaped oscillating tunnel was used to obtain detailed velocity profile under oscillatory boundary layer. In addition to these two experimental results, k-ω turbulence modeling is conducted for tsunami wave boundary layer as listed in Table 1. These are hypothetical shoaling tsunami cases investigated by Tinh and Tanaka [7]. Figure 4 indicates conditions of the whole test cases summarized in Table 1.
Most of the data have sufficiently high Reynolds number to fall in the turbulent regime. In addition, the present computational cases using k-ω model range from moderate to sufficiently high values of Re and a m /k s , due to an input of longer period of tsunami wave as seen in Table 1. In the numerical computation of k-ω model, the boundary layer equation is solved by finite difference method using an implicit scheme [31]. Sufficiently fine mesh size is employed in the immediate vicinity of the bottom, with gradually increasing mesh spacing in the vertical upward direction [32].
Considering the dimensionless quantity defined in Equation (15), the whole data sets are plotted in Figures 5 and 6 using the common vertical coordinate ζ.
It is clearly observed in Figures 5 and 6 that all data sets indicate overshooting velocity profile at ζ = 0.331, below which the velocity reduces gradually following logarithmic law with different velocity gradient. In addition, it is evident that Jensen [30] and the present numerical simulations are more suitable for detecting the exact elevation of overshooting height due to their finer interval of measurements and simulations.
In order to validate Equation (15), velocity profile data under wave crest ( = 0 deg) were collected as summarized in Table 1. In both laboratory experiments by Sleath [8] and Jensen [30], a U-shaped oscillating tunnel was used to obtain detailed velocity profile under oscillatory boundary layer. In addition to these two experimental results, k-turbulence modeling is conducted for tsunami wave boundary layer as listed in Table 1. These are hypothetical shoaling tsunami cases investigated by Tinh and Tanaka [7]. Figure 4 indicates conditions of the whole test cases summarized in Table 1. Most of the data have sufficiently high Reynolds number to fall in the turbulent regime. In addition, the present computational cases using k-model range from moderate to sufficiently high values of and / , due to an input of longer period of tsunami wave as seen in Table 1.     In the numerical computation of k-model, the boundary layer equation is solved by finite difference method using an implicit scheme [31]. Sufficiently fine mesh size is employed in the immediate vicinity of the bottom, with gradually increasing mesh spacing in the vertical upward direction [32].
Considering the dimensionless quantity defined in Equation (15), the whole data sets are plotted in Figures 5 and 6 using the common vertical coordinate .
It is clearly observed in Figures 5 and 6 that all data sets indicate overshooting velocity profile at  = 0.331, below which the velocity reduces gradually following logarithmic law with different velocity gradient. In addition, it is evident that Jensen [30] and the present numerical simulations are more suitable for detecting the exact elevation of overshooting height due to their finer interval of measurements and simulations.   For further investigation on wave boundary layer thickness, the height of the overshooting from these data sources is plotted in Figure 7. The variable on the abscissa is the roughness Reynolds number which can be utilized for classifying the flow regime under turbulent conditions into smooth and rough as shown previously [20].   For further investigation on wave boundary layer thickness, the height of the overshooting from these data sources is plotted in Figure 7. The variable on the abscissa is the roughness Reynolds number which can be utilized for classifying the flow regime under turbulent conditions into smooth and rough as shown previously [20].  Figure 7. Wave boundary layer thickness in terms of Equation (15).
The data sets listed in Table 1 fully encompass smooth, transitional, and rough regimes ( Figure  7). In addition, the constant defined in Equation (15) can be determined to be 0.331. Hence, we obtain the following equation, which is half theoretical and half empirical. * = 0.331 (22) It is important to note that in the theoretical treatment to derive Equation (22), the flow regime is assumed to be simply turbulent, and thus the classification into smooth and rough is not required. Hence it is concluded that Equation (22) is valid for both smooth and rough turbulent regimes. This is a distinct difference as compared with existing formulas from Equations (8)-(13).

Proposal of New Calculation Formulas for Each Flow Regime
Although Equation (22) has superiority over the existing expressions in the sense that it is valid for a turbulent regime regardless of bottom roughness condition, its convenience in practical applications might be improved by deriving an equation for each flow regime expressed in terms of Re and / . Therefore, friction coefficient under wave motion can be substituted into Equation (22). Here the following explicit form of wave friction coefficient for smooth bottom [20] and rough bottom [33] proposed by the one of the authors is used. The data sets listed in Table 1 fully encompass smooth, transitional, and rough regimes (Figure 7). In addition, the constant defined in Equation (15) can be determined to be 0.331. Hence, we obtain the following equation, which is half theoretical and half empirical.
It is important to note that in the theoretical treatment to derive Equation (22), the flow regime is assumed to be simply turbulent, and thus the classification into smooth and rough is not required. Hence it is concluded that Equation (22) is valid for both smooth and rough turbulent regimes. This is a distinct difference as compared with existing formulas from Equations (8)-(13).

Proposal of New Calculation Formulas for Each Flow Regime
Although Equation (22) has superiority over the existing expressions in the sense that it is valid for a turbulent regime regardless of bottom roughness condition, its convenience in practical applications might be improved by deriving an equation for each flow regime expressed in terms of Re and a m /k s . Therefore, friction coefficient under wave motion can be substituted into Equation (22). Here the following explicit form of wave friction coefficient for smooth bottom [20] and rough bottom [33] proposed by the one of the authors is used.
By substituting these two friction factors into Equation (22), a convenient formula for wave boundary thickness can be obtained for smooth and rough turbulent flows separately.
The equations obtained in the present study are illustrated in Figures 2 and 3. In the range of higher Re and a m /k s corresponding to tsunami-scale wave motion, the present equations show good agreement with Williams and Fuhrman [5].

Full-Range Equation
Sana and Tanaka [11] proposed a full-range equation applicable to all flow regimes spanning laminar, smooth turbulent, and rough turbulent. However, it is shown in this study that the equations for smooth turbulent and rough turbulent proposed by Sana and Tanaka [11] need to be replaced by Equations (25) and (26), respectively, to be more accurate up to higher values of Re and a m /k s . Accordingly, computation method of δ for an arbitrary flow regime including transitional regime is as follows [11]. δ where f 1 and f 2 are given by in which R 1 is defined by In order to apply Equation (27), the expressions for each flow regime, Equations (5), (25), and (26), will be substituted, rather than Equations (9) and (12). The result obtained by the full-range equation is illustrated in Figure 8. It is seen that highly smooth transitional behavior can be obtained from Equation (27). However, as compared with the previous diagram proposed by Sana and Tanaka [11], the transition from laminar to smooth turbulent flow is more abrupt.
In recent years, experimental studies have been carried out for boundary layers in an acceleration-skewed oscillatory flow (e.g., Suntoyo et al. [17], van der A et al. [34], and Yuan and Wang [35]). The equations proposed herein for a sinusoidal motion can be applied to acceleration-skewed oscillatory waves as well by introducing an equivalent sinusoidal orbital amplitude, a c , for a m as proposed by Suntoyo et al. [17] and van der A et al. [34]. In recent years, experimental studies have been carried out for boundary layers in an acceleration-skewed oscillatory flow (e.g., Suntoyo et al. [17], van der A et al. [34], and Yuan and Wang [35]). The equations proposed herein for a sinusoidal motion can be applied to accelerationskewed oscillatory waves as well by introducing an equivalent sinusoidal orbital amplitude, ac, for am as proposed by Suntoyo et al. [17] and van der A et al. [34].

Application of the New Formula to Demarcation of Friction Law under Tsunami
Tinh and Tanaka [7] proposed equations and a diagram demarcating wave friction zone and steady friction zone under shoaling tsunami. Using Equation (12), Tinh and Tanaka [7] obtained the following criteria for demarcating friction laws by assuming ℎ = . ℎ = 0.111 ( ) 0.754 (31) Meanwhile, using the newly obtained equation for boundary layer thickness for rough turbulent regime, Equation (26), an alternative equation is obtained by assuming ℎ = .
By equating wave friction factor and steady current friction factor , Tinh and Tanaka [7] proposed a demarcation in terms of friction factor.
In addition, Tanaka et al. [36] empirically determined the following criterion for frictional factor demarcation under tsunami.
Equations (31)- (34) are plotted in Figure 9 for comparison. It is observed in Figure 9 that Equation (31) is very far from other criteria. However, by replacing the previous boundary layer thickness formula with the improved one, newly obtained criterion in terms of boundary layer

Application of the New Formula to Demarcation of Friction Law under Tsunami
Tinh and Tanaka [7] proposed equations and a diagram demarcating wave friction zone and steady friction zone under shoaling tsunami. Using Equation (12) Equations (31)- (34) are plotted in Figure 9 for comparison. It is observed in Figure 9 that Equation (31) is very far from other criteria. However, by replacing the previous boundary layer thickness formula with the improved one, newly obtained criterion in terms of boundary layer thickness (Equation (32)) gives a result which is very close to Equations (33) and (34), especially at a higher value of a m /k s . It is interesting to note that the field data of Lacy et al. [6] during the 2010 Chilean tsunami can successfully be classified into the wave friction zone. (32)) gives a result which is very close to Equations (33) and (34), especially at a higher value of / . It is interesting to note that the field data of Lacy et al. [6] during the 2010 Chilean tsunami can successfully be classified into the wave friction zone. Since the equations proposed herein are based on linear wave theory, modification is required when we apply them to real tsunami events, in which wave motion is usually transient. Larsen and Fuhrman [37,38] proposed a method to use time-varying boundary layer thickness by replacing the orbital amplitude with the cumulative distance traveled by a free-stream particle.

Conclusions
In this study, investigation on wave boundary layer thickness was carried out, aiming to increase the applicability of the existing formulas to tsunami waves. According to comparisons among three existing formulas, it was found that Sana and Tanaka's [11] equation gives different behavior as compared with others at higher Re and / . Therefore, a new calculation formula was proposed based on the theoretical consideration using an eddy viscosity coefficient. By applying newly obtained equations for smooth and rough turbulent regimes, a full-range equation was successfully reconstructed, which enabled us to calculate wave boundary layer thickness automatically without judgment of flow regime for a set of given input parameters. In addition, the new equation can be successfully applied for demarcating friction factor under tsunami. The findings of this study will improve the prediction of sediment movement under long waves by providing more accurate boundary layer thickness under variable field conditions.
Author Contributions: H.T. wrote most of the manuscript and contributed to the interpretation of the results. N.X.T. and A.S. performed the experimental data analysis and numerical simulation as well as prepared for the figures. All authors have read, discussed, and agreed to the published version of the manuscript.
Funding: This research was funded by the Taisei Foundation.

Acknowledgments:
The authors would like to express their sincere gratitude to the financial support from the Taisei Foundation.  Since the equations proposed herein are based on linear wave theory, modification is required when we apply them to real tsunami events, in which wave motion is usually transient. Larsen and Fuhrman [37,38] proposed a method to use time-varying boundary layer thickness by replacing the orbital amplitude with the cumulative distance traveled by a free-stream particle.

Conclusions
In this study, investigation on wave boundary layer thickness was carried out, aiming to increase the applicability of the existing formulas to tsunami waves. According to comparisons among three existing formulas, it was found that Sana and Tanaka's [11] equation gives different behavior as compared with others at higher Re and a m /k s . Therefore, a new calculation formula was proposed based on the theoretical consideration using an eddy viscosity coefficient. By applying newly obtained equations for smooth and rough turbulent regimes, a full-range equation was successfully reconstructed, which enabled us to calculate wave boundary layer thickness automatically without judgment of flow regime for a set of given input parameters. In addition, the new equation can be successfully applied for demarcating friction factor under tsunami. The findings of this study will improve the prediction of sediment movement under long waves by providing more accurate boundary layer thickness under variable field conditions.